- Email: [email protected]
Data-driven thermally-induced error compensation method of high-speed and precision five-axis machine tools
Mechanical Systems and Signal Processing 138 (2020) 106538
Contents lists available at ScienceDirect
Mechanical Systems and Signal Processing journal homepage: www.elsevier.com/locate/ymssp
Data-driven thermally-induced error compensation method of high-speed and precision five-axis machine tools Jialan Liu a,b, Chi Ma a,b,⇑, Shilong Wang a,b a b
College of Mechanical Engineering, Chongqing University, Chongqing 400044, China State Key Laboratory of Mechanical Transmission, Chongqing University, Chongqing 400044, China
a r t i c l e
i n f o
Article history: Received 7 September 2019 Received in revised form 8 November 2019 Accepted 22 November 2019
Keywords: Five-axis machine tool Thermally-induced error Homogeneous transformation Error compensation
a b s t r a c t The data-driven thermal error compensation method of high-speed and precision five-axis machine tools was proposed based on the homogeneous transformation. The compensation component was obtained by analyzing the error transmission chain of machine tools and was expressed as the transmission matrix consisting of thermal error terms of linear axes and spindle system according to the differential movement of the compensation axis. From the view of thermal error generation mechanism, the thermal error of the linear axis was formulated as the product of the polynomial function with time as its independent variable and the polynomial function with position as its independent variable. Then the thermal errors of the linear axes were decomposed and identified from the measured positioning error. The thermal error of the spindle system was identified by the thermal characteristic experiments. The compensation component in each direction is calculated by substituting the identified thermal error terms of the linear axes and spindle system into the data-driven thermal error model. Finally, to demonstrate the effectiveness of the proposed method, the thermal error at a new working condition was predicted, and then the error compensation and actual machining were carried out. The results show that the machining error is reduced by more than 85% and 37% with the present error compensation compared with that without compensation and that without traditional error compensation, respectively. This research sheds new light on both the generation mechanism and the compensation method of the thermally-induced error of five-axis machine tools. Ó 2019 Elsevier Ltd. All rights reserved.
1. Introduction The demand of modern manufacturing to five-axis CNC machine tools is high-precision and high-speed machining, resulting from the pursuit of the machining quality and efficiency. High-speed and precision five-axis CNC machines tools are used to machine complex and precision parts in the automotive, aerospace, and military industry because it enables complex spatial movement between the tool and the workpiece. In particular, the use of high-speed and precision five-axis machine tools improves the machining accuracy and efficiency greatly. However, the high-speed and precision machine tools are under the combined effect of the geometric, dynamic, and thermal errors. Moreover, the thermally-induced error of the precision machine tools, caused by the uneven temperature distribution, accounts for up to 75% of the total volumetric error of machined workpieces [1]. So the precision of parts will be increased significantly if the thermally-induced error of ⇑ Corresponding author. E-mail addresses: [email protected] (J. Liu), [email protected] (C. Ma), [email protected] (S. Wang). https://doi.org/10.1016/j.ymssp.2019.106538 0888-3270/Ó 2019 Elsevier Ltd. All rights reserved.
2
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Nomenclature O1 O2 O3 O4 x y z b c s t B C S X Y Z
the intersection of B- and C-axes origin of the coordinate system b origin of the coordinate system c origin of the coordinate system s Local coordinate system fixed to X-axis Local coordinate system fixed to Y-axis Local coordinate system fixed to Z-axis Local coordinate system fixed to B-axis Local coordinate system fixed to C-axis Local coordinate system fixed to spindle system Coordinate system of tool nose B-axis C-axis S-axis X-axis Y-axis Z-axis a rotational angle of C-axis b rotational angle of B-axis kx angular errors in the compensation movement process ky angular errors in the compensation movement process kz angular errors in the compensation movement process rx translational errors in the compensation movement process ry translational errors in the compensation movement process rz translational errors in the compensation movement process dw change in coordinate system caused by differential transformation t T M2 friction torque between nut and screw shaft rotational speed of the screw shaft n2 Ph lead of screw shaft; V feed rate b0 helix angle of raceway Me friction torque caused by the external load Mg friction torque caused by geometric slip of balls mb and ma coefficients related to eccentricity of contact deformation ellipse # parameter related to material properties P q sum principle curvatures M gy gyroscopic moment F cq centrifugal force, F iq and F eq shear loads Q iq and Q eq normal loads A1 and A2 x and z coordinates of inner and outer groove curvature centers M1 friction torque of the bearings friction torque caused by the external loads Ml a the contact angle of the bearing. Max the greater value between 0:9F a =tana 0:1F r and F r Fa and Fr bearing’s axial and radial loads vo lubricant kinematic viscosity at operating temperature N rotational speed of bearing fo constant related to the bearing model and lubrication type Re Reynolds number Pr Prandtl number w angular velocity v viscosity of air diameter of the screw shaft d0 n1 rotational speed of screw shaft x and t the position and the time L length of screw shaft TL temperature of screw shaft C 1 and C 2 Coefficients related to temperatures
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
T ð xÞ position-dependent of temperature T ðt Þ time-dependent components of temperature T r ð xÞ position-component of rear bearing’s temperature T s ð0Þ initial temperature of nut’s temperature time-component of nut’s temperature T s ðt Þ oðxÞ geometric error Eðx; t; T Þ thermal error b0 , b1 , b2 , bm coefficients to be determined y output variable DY 1 and DY 3 translational errors in Y-direction cx thermal pitch / rotational angle of spindle system a0 , b0 , and c0 coefficients at x ¼ 0 gx and gy translational errors T f 1 and T f 2 temperatures of bearing #1 and bearing #2, specific heat capacity and mass of bearing cb and mb Ft preload r stress DL2 induced elongation DL3 pre-stretch capacity a0 thermal expansion coefficient y moving distance of Y-axis d1 offset of O2 in the x-direction of the reference system O1 d2 offset of O3 in the y-direction of the reference system O1 eYX squareness error between Y- and X-axes; Eðy; t; T Þ thermally-induced error of Y-axis dxy ; dyy ; dzy ; exy ; eyy and ezy six position-dependent geometric errors of Y-axis x moving distance of X-axis dxx ; dyx ; dzx ; exx ; eyx and ezx six position-dependent geometric errors of X-axis Eðx; t; T Þ thermal error of X-axis z moving distance of Z-axis dxz ; dyz ; dzz ; exz ; eyz and ezz six position dependent geometric errors of Z-axis eZX squareness error between Z- and X- axes eZY squareness error between Z- and Y- axes Eðz; t; T Þ thermal error of Z-axis L length of cutter d3 offset of O4 in z- direction of the reference system O1 E identity matrix cx , cy , and cz angular errors in x-, y-, and z- directions gx , gy , and gz translational errors in x-, y-, and z- directions ð X; cx Þ transformation matrixes around three x-axis R transformation matrixes around three y-axis R Y; cy R Z; cz transformation matrixes around three z-axis w homogeneous transformation matrix of the compensation movement of coordinate systems t to w t Tc Dx differential compensation movement of X-axis Dy differential compensation movement of Y-axis Dz differential compensation movement of Z-axis Da differential compensation movement of C-axis Db differential compensation movement of B-axis B the inner and outer groove curvature constant; a0 the initial contact angle D the diameter of balls h the angular displacement of the inner ring wq position angle of qth ball Ri the radius of the groove curvature center contact deformations of balls and grooves da and dr K i and K o contact stiffness aeq and aiq contact angles Fa axial load radial load Fr Mg moment load
3
4
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
X 1 and X 2 initial values of x and z coordinates of balls n1 rotational speed of the screw shaft Mv friction torque caused by lubricant viscous friction f1 coefficient related to the bearing structure equivalent load applied on bearings Fb dm pitch diameter of bearings L feature size k thermal conductivity of fluid Nu Nusselt number T0 temperature of ambient; h convective coefficient q density and the respectively; c specific heat capacity d0 denote diameter of screw shaft T0 reference temperature k thermal conductivity T ð0Þ temperature of heat source DT ðx; t Þ temperature rise initial temperature of front bearing’s temperature T f ð0Þ T f ð xÞ position-component of front bearing’s temperature T f ðt Þ time-component of front bearing’s temperature T r ð0Þ initial temperature of rear bearing’s temperature time-component of rear bearing’s temperature T r ðt Þ T s ðxÞ position-component of nut’s temperature At , Bt , and C t coefficients related to time ax , bx , and cx coefficients related position e the residue x1 , x2 , , and xm input variables T 1 , T 2 , , and T m critical temperatures cy thermal yaw DX 2 and DX 4 translational errors in X-direction R distance between two sensors on the test rod DX 1 and DY 1 distances in X- and Y- directions from probe of temperature sensor S1 to test rod g Fitting ability/Predictive ability Tn temperature of moving nut cn and mn specific heat capacity and mass of nut Fa axial load A and P denote cross-section area and axial force, DL1 actual elongation DL thermal error with pre-stretch
five-axis CNC machines tools is reduced or avoided. So far, the thermally-induced error investigation becomes a hot topic of current research activities. The error avoidance method, the thermal characteristic optimization method, and the error compensation method are widely applied to eliminate the thermally-induced error. The error avoidance method attempts to eliminate or reduce possible error sources at the design and manufacturing stages [2]. The influence of error sources on the system is reduced with the increase in the design and manufacturing accuracy of machine tools’ components, and then the strict temperature control, the vibration isolation, the airflow disturbances, and the environmental state control are employed to eliminate or reduce the influence of external error sources of machine tools. Although it can reduce the original error, there still exist great limitations in the improvement of the manufacturing and installation accuracy of machine tools. It is shown that the cost of the error avoidance method increases exponentially when the machining accuracy requirements are higher than a certain level. Then the error avoidance method fails in improving the machining accuracy of machine tools costly. When the error avoidance method encounters difficulties, researchers switch their interest to the thermal and dynamic characteristic analysis and optimization of machine tools. The bearing thermal model [3–4], the thermal-structure coupling characteristics analysis [5–6], the thermal contact resistance (TCR) [7–9], the optimization of the fits for the bearing joints [10], the hydrodynamics analysis of machine tools [11], the contact stiffness of the tool/holder interface [12], and the precision loss modeling of machine tools [13] are involved in the thermal characteristic analysis of machine tools. Moreover, the thermal characteristic analysis is not only conducted for the whole machine tools, but also is conducted for critical functional components of machine tools [14–20]. The thermal characteristic analysis method is effective to analyze the temperature rise and the thermal deformation of machine tools. But accurate heat source intensity and thermal boundary conditions
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
5
are difficult to be obtained because of the nonlinear temperature rise and thermal displacements of machine tools, resulting in the modeling error of the thermal characteristic. Moreover, the results obtained by the thermal characteristic analysis are difficult to be used in the subsequent error compensation because the numerical modeling is time-consuming. The above factors limit the wide application of the thermal characteristic analysis and optimization methods. The thermally-induced error compensation is then becoming more and more popular in the past few decades. Moreover, it is revealed that the thermally-induced error is reduced significantly when the compensation method is implemented [21– 23]. The main process includes such steps as thermal error measurement, identification, modeling, and compensation. The thermal error measurement and identification are realized by the laser interferometer, etc. [24–25]. The thermal error model with excellent predictive performance and strong robustness is the core of the thermal error compensation according to previous studies [26–27]. Then various algorithms are used to establish the correlation between the thermal error and measured temperatures. The most widely-used algorithms in the thermal error modeling are the multivariate linear regression analysis (MLRA) [28], the time series (TS) [29], the genetic algorithm (GA) [30], and the artificial neural network (ANN) [31]. The error compensation of linear axes is also addressed to reduce the thermally-induced error [32–35]. However, the robustness of the above models needs further verification because of the dynamic change in the thermal behaviors, resulting from complex and volatile working conditions. The above models are the empirical error models with temperatures as the input and with the thermal error as the output, and the modeling method is less compatible with thermal error data without revealing the error generation mechanism. Moreover, for the five-axis machine tools, the error generation mechanism of has not been fully investigated to construct the thermal error model. Instead, the thermally-induced error is correlated with critical temperatures empirically, which splits the relationship between the error generation mechanism and the modeling method. Furthermore, the thermal error model without physical significance is difficult to be used in practical application, including the thermal errors based on GA and ANN. Besides, the combined effects of thermal errors of the linear and rotational axes are considered simultaneously, leading to the failure in the significant improvement of the machining accuracy of five-axis machine tools. Then the integrated thermal error modeling and compensation methods with strong robustness and excellent predictive performance are necessary to be investigated. To realize the integrated thermal error modeling and compensation, the homogeneous transformation (HT) method is used [36–41]. The thermal characteristic experiments are conducted offline, and then the error terms are identified. The error terms are substituted into the transformation matrix and these error terms generally are the static errors at thermal equilibrium state in the above models, failing in the dynamic modeling and compensation. Moreover, the mapping relation between the thermally-induced error and critical temperatures is constructed empirically because it is believed that the thermal error is the direct reflection of the temperature rise. But the thermal error generation mechanism of the linear and spindle axes is somewhat complex because of the interaction among such factors as the structure parameters, the working conditions, the ambient temperature, the material properties, the thermal load intensity and position, and the cooling effect, leading to the deterioration of the robustness. The relationship between the temperature rise, the thermal expansion, the thermally-induced error, and the resultant positioning error of linear axes has not been explored. So there still exists shortcomings in the HT method, and it is difficult to realize accurate and robust error compensation. The main goal of this study is to reduce the thermally-induced error of five-axis machine tools. The remarkable difference from previous studies is the integration of the thermal errors of linear and rotational axes with the aid of the homogeneous transformation, and the thermal errors of the linear axes are identified with the product of the polynomial function with temperature as its independent variable and the polynomial function with position as its independent variable from the error generation mechanism of linear axes. The advantage of the data-driven thermal error compensation method is its excellent effectiveness and strong robustness at different working conditions. The most significant contribution is the realization of the thermally-induced error compensation by the compensation model with excellent predictive performance and strong robustness. The main principle of the proposed method is expressing the compensation components as the matrix of error terms of the linear and spindle axes by analyzing the error transmission chain of five-axis machine tools. In Section 2, the thermal error model was established by calculating the homogenous transformation matrixes of the workpiece chain and tool chain of the five axis CNC machine tools, and then the compensation component was obtained for each axis. In Section 3, the thermal error terms of the linear and spindle axes were identified by the experiments. In Section 4, the thermal error compensation was implemented for a five-axis CNC machine tool to verify the validity of the proposed method. In Section 5, the effects of the heat source and pre-stretch of the lead screw on the thermal error are discussed. In Section 6, the conclusions obtained from previous analysis were summarized.
2. Thermally-induced error modeling method of five-axis machine tool 2.1. Thermally-induced error transfer model A standard five-axis CNC machine tool is made up of two rotational axes of B and C and three translational axes of X, Y, and Z, as depicted in Fig. 1. In the machining process, the cutter is installed on the spindle system in Z-direction, and the workpiece is installed on C-axis. The settings of the coordinate system are listed in Table 1 to assess the thermallyinduced error of the five-axis CNC machine tool. Three translational axes X, Y, and Z are perpendicular with each other mutually. The directions of the local coordinate systems are the same as that of the five-axis machine tool’s coordinate system. The
6
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 1. Configuration schematic view of five-axis CNC machine tools.
Table 1 Definition of coordinate systems. Coordinate system
Physical significance
x y z b
Local Local Local Local
c
Local coordinate system fixed to C-axis
s
Local coordinate system fixed to spindle system Coordinate system of tool nose Coordinate system of workpiece
t w
coordinate coordinate coordinate coordinate
system system system system
Directions fixed fixed fixed fixed
to to to to
X-axis Y-axis Z-axis B-axis
Remarks
Same as the direction of coordinate system of machine tools Coordinate origin O2 Coordinate origin O3 Coordinate origin O4
origins of the coordinate systems b, c, and s are O2 , O3 , and O4 , respectively. Moreover, the point O2 is on the B-axis, the point O3 is at the center point of the top surface of C-axis, and the point O4 is at the center point of the front end surface of the cutter. The reference coordinate system O1 is the intersection of B- and C- axes, and the intersection O1 is on the axis of S-axis. The transmission path of the thermally-induced error of five-axis machine tools is shown in Fig. 2. Moreover, the total transmission paths are divided into two transfer kinematic chains. One of the transmission paths of the kinematic chain is expressed successively as the machined workpiece!rotational axis C!rotational axis B!translational axis Y!translational axis X!machine tool bed. The other transmission path of the kinematic chains is expressed as the tool!the spindle system S!the translational axis Z!machine tool bed. The error will be transferred to each axis by the above two kinematic chains, resulting in the deterioration of the machining accuracy. The transformation matrixes between different coordinate systems are expressed as
7
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 2. Schematic view of kinematic chain.
2 b cT
cos a sin a 0 d2
6 sin a 6 ¼6 4 0
3
2
cos b
0
sin b
cos a
0
0
1
6 0 d1 7 1 0 7y 6 7 bT ¼ 6 4 sin b 0 cos b 0 5
0
0
1
0
0
0
d2
3
2
1 0
0
0
3
60 1 0 y7 07 7x 7 6 7 yT ¼ 6 7 40 0 1 05 05
0
0 0
1
0
ð1Þ
1
where a denotes the rotational angle of C-axis; b denotes the rotational angle of B-axis; y denotes the moving distance of Y-axis; d1 denotes the offset of O2 in the x-direction of the reference system O1 ; and d2 denotes the offset of O3 in the y-direction of the reference system O1 . Then the actual homogeneous transformation matrix is corrected with the geometric and thermal errors of the Y-axis.
2 x 0 yT
eYX
1
6e 6 YX ¼6 4 0
0
32
1 0
0
0
32
0 0
0 1
0
0
0
0 1
0 eyy dxy 1 0 0 6 0 1 0 Eðy; t; T Þ 7 exy dyy 7 7 76
1
exy
1
0
0
0
3
32
ezy
1
76 6 07 76 0 1 0 y 76 ezy 76 76 4 5 0 0 1 0 54 eyy 1 0
1
0
0
76 dzy 54 0 0
1
0 1
0
0
1
0
ð2Þ
7 5
where eYX denotes the squareness error between Y- and X-axes; Eðy; t; T Þ denotes the thermally-induced error of Y-axis; and dxy ; dyy ; dzy ; exy ; eyy and ezy denote six position-dependent geometric errors of the Y-axis. The actual homogeneous transformation matrix is corrected with the geometric and thermal errors of the X-axis.
2 O1 0 x T
1 0 0
x
32
1
6 0 1 0 0 76 e 76 zx 6 ¼6 76 4 0 0 1 0 54 eyx 0
0
0 1
eyx dxx 1 0 0 Eðx; t; T Þ 7 6 0 exx dyx 7 7 76 0 1 0
1
0
3
32
ezx
exx
1
76 dzx 54 0
0
0
1
0
0 1
0
0
1
0
ð3Þ
7 5
where x denotes the moving distance of the X-axis; dxx ; dyx ; dzx ; exx ; eyx and ezx denote 6 position-dependent geometric errors of the X-axis; and Eðx; t; T Þ denotes the thermal error of the X-axis. The actual homogeneous transformation matrix is corrected with the geometric and thermal errors of the Z-axis.
2
1 6 0 6 O1 0 z T ¼ 6 4 eZX 0
0 1
eZY 0
32
32
3
32
1 ezz eyz dxz 1 0 0 0 eZX 0 1 0 0 0 76 7 6 6 0 eZY 0 7 1 exx dyz 7 7 76 0 1 0 0 76 ezz 76 0 1 0 76 7 76 76 1 0 54 0 0 1 z 54 eyz exz 1 dzz 54 0 0 1 Eðz; t; T Þ 5 0
1
0
0
0 1
0
0
0
1
0
0
0
ð4Þ
1
where z denotes the moving distance of Z-axis; dxz ; dyz ; dzz ; exz ; eyz and ezz denote 6 position dependent geometric errors of the Z-axis; eZX and eZY denote the squareness error between Z- and X- axes and that between the Z- and Y- axes, respectively; and Eðz; t; T Þ denotes the thermal error of the Z-axis.
8
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
According to the transfer relationship of the thermal and geometric errors between two kinematic chains, the following relationship is obtained. O1 0 z T
z sT
s tT
¼ Ox 1 T 0
x 0 yT
y bT
b cT
ð5Þ
c wT
In the ideal case, the transformation matrix of coordinate systems w to c is expressed as c wT
¼ ðOz 1 T 0
where ts T , sz T , and
0 z O1 T
2 t sT
¼
O1 0 x T
s tT Þ
z sT
x 0 yT
y bT
b cT
¼
t sT
s zT
0 z O1 T
O1 0 x T
denote the inverse matrix of the st T , zs T , and
1
0
0
0
O1 0 z T ,
x 0 yT
y bT
respectively. Then st T , sz T , and
0
0 z O1 T
are expressed as
3
7 s 1 6 60 1 0 07 ¼6 7 tT 40 0 1 L5 0
ð6Þ
b cT
ð7Þ
0 1
where L denotes the length of the cutter.
3 1 0 0 0 60 1 0 0 7 7 6 ¼6 7 4 0 0 1 d3 5 2
s zT
¼
z 1 sT
0
0
0
ð8Þ
1
where d3 denote the offset of O4 in z-direction of the reference system O1 . 0 z O1 T
¼ ðOz 1 T 0 Þ
1
ð9Þ
The transformation matrix of coordinate systems t to w is expressed as w t T
O1 0 O1 0 z s w c b y 0x ¼E¼w O1 T t T ¼ c T b T y T x T O1 T z T s T t T
ð10Þ
where E denotes the identity matrix. 2.2. Thermally-Induced error modeling principles The thermal error terms of the spindle system include the angular errors of cx , cy , and cz and the translational errors of gx , gy , and gz , as shown in Fig. 3. Then the transformation matrixes of the spindle system around three transitional axes are obtained.
Fig. 3. Schematic graph of thermal error components of spindle system.
9
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
2
1
0
0
6 0 cosc 6 x Rð X; cx Þ ¼ 6 4 0 sincx 0
coscx 0
0 coscz
6 sinc 6 z RðZ; cz Þ ¼ 6 4 0
ð11Þ
1
0
sincy
1
0
0
3
0 coscy
07 7 7 05
0
1
0
sincz
0
0
ð12Þ
3
0
07 7 7 1 05
0
0 1
coscz
0
3
07 7 7 05
sincx
0
2 coscy 6 0 6 R Y; cy ¼ 6 4 sincy 2
0
0
ð13Þ
The spindle system rotates continuously around three coordinate axes of a fixed coordinate system for three times, and then
2
coscy coscz
6 sinc sinc cosc þ cosc sinc 6 x z z x z R ¼ Rð X; cx ÞR Y; cy RðZ; cz Þ ¼ 6 4 coscx sincy coscz þ sincx coscx 0
coscy sincz
sincy
sincx sincy sincz þ coscx coscz
sincx coscy
sincz sincz þ sincx coscz
coscx coscy
0
0
0
3
07 7 7 05 1 ð14Þ
Based on the assumption of small errors, then coscx ¼ coscy ¼ coscz ¼ 1, sincx ¼ cx , sincy ¼ cy and sincz ¼ cz . If the terms higher than the second-order are ignored, then Eq. (14) is simplified as
2
1
6 c 6 R ¼6 z 4 cy 0
cz 1
3
cy 0 cx 0 7 7
cx
1
7 05
0
0
1
ð15Þ
Combining Eq. (15) with the translational error of the spindle system, the transformation matrix movement of the tool nose coordinate systems t to w is expressed as
2 w t Te
1
6 c 6 ¼6 z 4 cy 0
cz 1
cx 0
twT e of the error
3
cy gx gx gy 7 7 7 1 gz 5 0
ð16Þ
1
Then the error compensation is achieved by the differential motion of each axis, and the error generated in the compensation movement is ignored. The compensation movement is the inverse movement of the error movement. The coincidence of the cutter tool nose coordinate system t and the workpiece coordinates system w is realized when the compensation movement is carried out. Namely, w w t Tct Te
¼E
ð17Þ
where twT c denotes the homogeneous transformation matrix of the compensation movement of the coordinate systems t to w. By using the differential transformation principle, the total differential relationship between coordinate systems t and w is obtained.
dw t T ¼
@w @wT @wT @wT @wT t T Dx þ t Dy þ t Dz þ t Da þ t Db @x @y @z @a @b
ð18Þ
where Dx, Dy, Dz, Da and Db denote the differential movement of the compensation of the translational axes of X, Y, and Z and the rotational axes of C and B. Moreover, the partial differential is expressed as y @w c b y 0 dx T O1 0 z s t T ¼w T T T c T bT y T x T @x dx z s t
ð19-aÞ
10
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538 y
dO1 T @w c b t T ¼w c T bT y T @y dx
0
0 O1 z s x O1 T z T s T t T
ð19-bÞ
@w dOz 1 T z s c b x y t T ¼w T T c T b T x T y T O1 T @z dz s t
ð19-cÞ
dcb T b y x O1 z s @w t T ¼w T T T T T T c T @a da y x O1 z s t
ð19-cÞ
dby T y x O z s @w c t T ¼w T T 1T T T c T bT @b db x O1 z s t
ð19-dÞ
The homogeneous transformation matrix of the compensation movement between the coordinate systems t and w is expressed as
2
1 6 k z 6 w t T c ¼6 4 ky 0
kz
ky
1 kx
kx 1
0
0
rx 3 ry 7 7 7 rz 5
ð20Þ
1
where kx , ky , and kz denote the angular errors in the compensation movement process, and rx , ry , and rz denote the transw lational errors in the compensation movement process. It is revealed that w t T c has an identical form with t T e . According to the error compensation principle, the following differential relationship is obtained. w t Tc
w ¼w t T þ dt T
ð21Þ
where dw t T denotes the change in the coordinate system caused by the differential transformation. Then the following relationship is obtained according to Eq. (21).
kx ¼ cx ; ky ¼ cy ; kz ¼ cz
ð22-aÞ
rx ¼ gx ; ry ¼ gy ; rz ¼ gz
ð22-bÞ
It is shown that the unique solution cannot be obtained by Eq. (22). The angular error cz of the tool cutter has no effect on machining accuracy of CNC machine tools because the rotational movement of the tool cutter relative to the spindle coordinate system z is symmetrical. Therefore, it is not necessary to compensate the thermal error in Z-direction of the workpiece coordinate system w and the cutter tool nose coordinate system t. Then the thermal error compensation components are obtained as
2 @ w T 31 t 3 2 3 2 3 gx rx Dx 7 6 @x @w t T 7 7 6 g 7 6 6 6r 7 7 6 @y 6 Dy 7 6 6 y7 y7 7 6 7 6 @w T 7 7 6 7 6 6 Dz 7 ¼ 6 gz 7 ¼ 6 6r 7 t 7 7 6 7 6 6 @z 7 6 z 7 6 7 6 7 6 6 7 w 7 4 Da 5 4 cx 5 6 6 @@t aT 7 4 kx 5 5 4 cy ky Db @w T 2
ð23Þ
t
@b
It is seen that the unique solution is obtained by Eq. (23). Then the relationship between the error compensation components Dx, Dy, Dz, Da, and Db are obtained, and then the thermal errors gx , gy , gz , cx , and cy of the spindle system are then obtained. Moreover, it is noticed that
@w T @w T @w T @w T t t t , @y , @z , @t a , @x
and
@w T t @beta
are closely related to geometric errors and thermal errors
of the linear axes. Then the following error terms are to be identified, including the geometric error terms dxx ; dyx ; dzx ; exx ; eyx , ezx and the thermal error term Eðx; t; T Þ of X-axis, the geometric error terms dxy ; dyy ; dzy ; exy ; eyy , ezy , and eYX , and the thermal error term Eðy; t; T Þ of Y-axis, the geometric error terms dxz ; dyz ; dzz ; exz ; eyz , ezz , eZX , eZY and the thermal error term Eðz; t; T Þ of Z-axis, and the thermal error terms cx , cy , cz , gx , gy and gz of the spindle system. The elimination of thermal error terms of the linear and spindle axes is the focus, including Eðx; t; T Þ, Eðy; t; T Þ, Eðz; t; T Þ, cx , cy , cz , gx , gy and gz . Then the identification of the geometric errors of exx ; eyx ; ezx , exy ; eyy ; ezy , exz ; eyz ; ezz , eYX , eZX , and eZY will be not studied anymore in the following section, and the identification of the thermal error terms of the linear and spindle axes, including Eðx; t; T Þ, Eðy; t; T Þ, Eðz; t; T Þ, cx , cy , cz , gx , gy , and gz will be investigated. As is shown, the compensation components in different directions are closely related to the position, time, temperature, the offsets d1 , d2 , and d3 , and the length L of the cutter. Moreover, the compensation component of each axis can be obtained according to Eq. (23) if the thermal error terms of the linear and spindle axes are identified.
11
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
2.3. Thermally-Induced error modeling method 2.3.1. Thermally-Induced error modeling of linear axis The thermal error modeling method of the linear axis will be studied from the view of the thermal error generation mechanism. The temperature field of the screw shaft is analyzed, and then the thermal error model of the linear axis will be proposed based on the thermal field modeling. (1) Heat Generation of Screw Nut Pair The heat generated by the motion of the screw nut pair is expressed as
Q2 ¼
2p n2 M2 60
ð24Þ
where M 2 denotes the friction torque between the nut and screw shaft; and n2 denotes the rotational speed of the screw shaft and can be expressed as
n2 ¼
V Ph
ð25Þ
where P h denotes the lead of the screw shaft; and V denotes the feed rate. The friction torque of the screw shaft pair is expressed as
M ¼ 2zðMe þ Mg Þcosb0
ð26Þ
0
where b denotes the helix angle of the raceway; and Me and Mg denote the friction torques caused by the external load and the geometric slip of balls, respectively. According to Ref. [20], Me and Mg are expressed as
sffiffiffiffiffiffiffiffiffiffiffiffiffi 4Q 4iq 3 P Me ¼ mb # q
ð27aÞ
vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi u 16Q 5iq 3 f m2a u t Mg ¼ 0:08 P R ð# qÞ2
ð27bÞ
where ma and mb denote the coefficients related to the eccentricity of the contact deformation ellipse, respectively; f denotes the coefficient of sliding friction; # denotes the parameter related to the material properties and is expressed as P q denotes the # ¼ 28 2 , wherein l and E denote the Poisson’s ratio and modulus of elasticity, respectively; and 3
1l 1l 1þ 2 E1 E2
sum principle curvatures of the contacting bodies at the contact point on the orthogonal plane and can be expressed as
X
q ¼ q11 þ q12 þ q21 þ q22
ð28Þ
a For the groove surface of the screw side, q11 ¼ q12 ¼ D2 , q21 ¼ 1R and q22 ¼ 2cosbcos , in which R denotes the raceway radius. dDcosb a For the groove surface of the nut side, q11 ¼ q12 ¼ D2 , q21 ¼ 1R and q22 ¼ 2cosbcos . dDcosb
The load imposed on the balls is shown in Fig. 4(a). As is shown, the ball is under the combined effect of the gyroscopic moment M gy , centrifugal force F cq , the shear loads F iq and F eq and the normal loads Q iq and Q eq . When the rotational speed of
Fig. 4. Quasi-static mechanics analysis of ball screw nut.
12
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
screw shaft is relatively high, the positions of the balls and the curvature centers of inner and outer groove are shown in Fig. 4(b). A1 and A2 denote the x and z coordinates of the inner and outer groove curvature centers and can be obtained according to Ref. [5].
(
A1 ¼ BDsina0 þ da þ hRi coswq
ð29Þ
A2 ¼ BDcosa0 þ dr coswq
where B denotes the inner and outer groove curvature constant; a0 denotes the initial contact angle; D denotes the diameter of the balls; h denotes the angular displacement of the inner ring; wq is the position angle of the qth ball; Ri denotes the radius of the groove curvature center; and da and dr denote contact deformations of the balls and grooves and are obtained according to the Pythagorean theorem and the geometric relationship shown in Fig. 4(b).
1=2 da ¼ V 21 þ V 22 ðf i 0:5ÞD h i1=2 ðf o 0:5ÞD dr ¼ ðA1 V 1 Þ2 þ ðA2 V 2 Þ2
ð30Þ
where f i and f o denote the curvature radius coefficients of inner and outer grooves, respectively; and V 1 and V 2 denote the distance variation between the curvature center of the inner ring and balls, respectively. The relationship between the center of the qth ball and the curvature centers of inner and outer grooves is expressed as
(
ðA1 X 1 Þ2 þ ðA2 X 2 Þ2 ¼ ½ðf i 0:5ÞD þ da X 21
þ
X 22
¼ ½ðf o 0:5ÞD þ dr
2
2
ð31Þ
where X 1 and X 2 denote the initial values of the x and z coordinates of balls, respectively. The qth ball is under the combined effect of the gyroscopic moment M gy , centrifugal force F cq , the shear loads F iq and F eq and the normal loads Q iq and Q eq . Then, the force balance equations of the qth ball is expressed as M gq cos eq d M gq sin eq iq d
Q oq sinaeq Q iq sinaiq Q oj cosaeq Q iq cosa
a cosaiq ¼ 0 a sinaiq þ F cq ¼ 0
ð32Þ
3=2 where Q iq ¼ K i d3=2 iq ; Q oq ¼ K o doq , wherein K i and K o denote the contact stiffness of the ball / nut raceway and ball / screw
raceway; aeq and aiq denote the contact angles of the balls / nut groove and balls / screw shaft groove, respectively. The principal curvature function is expressed as
F ðqÞ ¼
j2 þ 1 R 2C ðj2 þ 1ÞR
ð33Þ
P 0:636 P P P q q q where j ¼ 1:0339 P q1 ; C ¼ 1:5277 þ 0:6023ln P q1 ; and R ¼ 1:0003 þ 0:5968ln P q1 , wherein q1 ¼ q11 þ q22 2 2 2 P and q2 ¼ q12 þ q21 . The axial load F a , radial load F r , and moment load M g of the inner ring are expressed as
8 z P > > Fa Q iq sinaiq ¼ 0 > > > q¼1 > > > < z P Q iq cosaiq ¼ 0 Fr > q¼1 > > > > z P > > > Q iq Ri sinaiq ¼ 0 : Mg
ð34Þ
q¼1
where z denotes the number of balls. The model of the screw nut pair is W4009SA-6D-C5Z12 and is manufactured by the NSK Company. The basic parameters of the screw nut pair are as follows: the total length is 1269 mm, the diameter is 40 mm, the length of the nut is 225 mm, the maximum stroke is 700 mm and the lead is 12 mm. The rotational speeds of the screw shaft are 416.7 r/min, 666.7 r/min, 1000 r/min and 1250 r/min when the feed rates are 5 m/min, 8 m/min, 12 m/min and 15 m/min, respectively. The heat generation of the screw nut pair under different feed rates is obtained, seen in Fig. 5. As is shown, the heat generation of the screw nut pair increases with the feed rate. (2) Heat Generation of Bearings According to Ref. [6], the heat generated by the bearings is expressed as
(
Q 1 ¼ 1:047 104 n1 M1 M1 ¼ Ml þ Mv
ð35Þ
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
13
Fig. 5. Heat generation of screw nut pair.
where M 1 denotes the friction torque of the bearings; n1 denotes the rotational speed of the screw shaft; and M l and M v denote the friction torques caused by the external loads and the lubricant viscous friction, respectively. The friction torque M l caused by the external load is expressed as
M l ¼ f 1 F b dm
ð36-aÞ
F b ¼ maxð0:9F a =tana 0:1F r ; F r Þ
ð36-bÞ
where f 1 denotes the coefficient related to the bearing structure; F b denotes the equivalent load applied on bearings; dm denotes the pitch diameter of bearings; Fa and Fr denote bearing’s axial and radial loads, respectively; max denotes the greater value between 0:9F a =tana 0:1F r and F r ; and a denotes the contact angle of the bearing. The friction torque M v caused by lubricant viscous friction is expressed as
(
Mv ¼ 107 f o ðv o nÞ2=3 dm Mv ¼ 160 10
7
f o dm
3
3
v o n P 2000 v o n < 2000
ð37Þ
where v o denotes the lubricant kinematic viscosity at operating temperature; n denotes the rotational speed of the bearing; and f o denotes the constant related to the bearing model and lubrication type. f o is taken as 1 for single oil-mist or oil-gas lubrication. f o is taken as 2 for oil lubrication and is taken as 4 for grease lubrication of double column. The quasi-static analysis method, which is similar with the method shown in Ref. [3] to Ref. [6], is used to compute the heat generation and the friction torques of the bearings and the ball screw. The nonlinear balance equation set of bearings, including the balance equations of the geometrical, physical, rolling element, and bearing rings, are established. Then the Newton-Rapson algorithm is used to solve nonlinear balance equation set and the heat generation and the friction torques of the bearings are obtained. The parameters of the bearings are as follows. The inner diameter, the outer diameter, the thickness, and the contact angle, and the curvature radius coefficient of the grooves are 40 mm, 90 mm, 20 mm, 60°, and 0.55, respectively. Then the heat generation of bearings is obtained, seen in Fig. 6. It is depicted that the heat generation increases with the feed rate and the preload. (3) Convective coefficient According to Ref. [42], the convective coefficient is expressed as
14
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 6. Heat generation of bearings.
h¼
Nu k L
ð38Þ
where L denotes the feature size; k denotes thermal conductivity of fluid; and Nu denotes the Nusselt number. The heat transfer mode of moving components, including the screw shaft, belongs to forced convection. The criterion equation for the forced convection is expressed as
Nu ¼ 0:133Re2=3 Pr1=3
ð39Þ
where Re and Pr denote the Reynolds number and Prandtl number, respectively. The Reynolds number Re is expressed as
Re ¼ wd0 =v 2
ð40Þ
where w, v , and d0 denote the angular velocity, the viscosity of air, and the diameter of the screw shaft, respectively. The angular velocity is expressed as
w ¼ 2pn2
ð41Þ
where n1 denotes the rotational speed of the screw shaft. The forced convective coefficients of the screw shaft under different feed rates are obtained according to the above analysis, as shown in Table 2. It is seen that the forced convective coefficient increases with the fee rate. (4) Thermal Field Modeling The cross-section area and the perimeter of the shaft along the longitudinal direction do not change by assuming that the screw shaft is an extension rod with a circular cross-section. The temperature of the screw shaft at the axial position x and time t is T ðx; tÞ. Then the temperature rise of the screw at the different axial positions is expressed as T ðx; tÞ T 0 . According to Fig. 7, the heat transfer equation is obtained as
Table 2 Forced convective coefficients of screw shaft. Feed rate (m/min)
Rotational speed (r/min)
Forced convective coefficients (W= m2 C )
5 8 12 15
416.7 666.7 1000 1250
70.1 75.6 81.2 86.4
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
15
Fig. 7. Heat transfer of the shaft.
@ 2 T ðx; tÞ qc @T ðx; tÞ 4h þ dx ¼ ðT ðx; t Þ T 0 Þ @x2 k @t kd0
ð42Þ
where T 0 denotes the temperature of the ambient; h denotes the convective coefficient; and a0 denotes the thermal expansion coefficient; q and c denote the density and the specific heat capacity, respectively; d0 and T 0 denote the diameter of the screw shaft and the reference temperature, respectively; k denotes the thermal conductivity; x and t denote the position and the time, respectively. When the linear axis reaches the thermal equilibrium state, Eq. (42) is written as 2
d T ð xÞ 4h ¼ ð T ð xÞ T 0 Þ dx2 ad0 Assuming that
4h
ad0
ð43Þ
¼ m2 , then Eq. (43) is expressed as
2
d T ð xÞ ¼ m2 T ðxÞ dx2
ð44Þ
If the thermal expansion coefficient a is a constant, Eq. (42) is a quadratic linear differential equation with constant coefficients. Then the solution of Eq. (44) is expressed as
T ðx; t Þ ¼ C 1 emx þ C 2 emx
ð45Þ
Moreover, the following relationships are satisfied by Eq. (45). (1) T ¼ T ð0Þ when x ¼ 0, wherein T ð0Þ denotes the temperature of the heat source. ¼ aL T L , wherein L and T L denote the length and the temperature of the screw shaft, respectively. But T L is (2) k dT dx x¼L unknown, then it is necessary to add another condition to determine C 1 and C 2 . Namely, (3) T ¼ T L when x ¼ L. Then the following relationship is obtained.
T ð0Þ ¼ C 1 þ C 2
ð46Þ
C 1 memL C 2 memL ¼ aT L
Then C 1 and C 2 are obtained according to Eq. (46), and then the temperature distribution in the axial-direction is expressed as
T ðx; t Þ ¼
ðkm þ hÞemðLxÞ þ ðkm hÞemðLxÞ T ð0Þ ðkm þ hÞemL þ ðkm hÞemL
ð47Þ
The temperature rise stage is the none-equilibrium state. The efficiency of the heat conduction is far higher than that of the convection. Then the heat transfer equation is simplified as
@ 2 T ðx; tÞ qc @T ðx; tÞ dx ¼ @x2 a @t
ð48Þ
By solving Eq. (48), the temperature response is then obtained.
"
2 T ðx; t Þ ¼ T ð0Þ 1 pffiffiffiffi
p
Z px ffiffiffiffi a 2
qct
#
e
u2
du
ð49Þ
0
Then the temperature model of the screw shaft under the effect of a single heat source is obtained according to Eq. (45) and Eq. (49).
16
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
T ðx; tÞ ¼ T ð0Þ
ðkm þ hÞemðLxÞ þ ðkm hÞemðLxÞ 2 1 pffiffiffiffi ðkm þ hÞemL þ ðkm hÞemL p
Z px ffiffiffiffi a 2
qct
! 2
eu du
þ T 0 ¼ T ð0Þ T ðxÞ T ðt Þ þ T 0
ð50Þ
0
where T ðxÞ and T ðt Þ denote the position-dependent and time-dependent components of the temperature, respectively. It is necessary to know the parameters in Eq. (50) to get the thermal field. Then the parameters in Table 3 are used to compute the thermal field at the thermal equilibrium state, seen in Fig. 8. It is shown that the temperature increases with time and that the temperature does not change significantly with position at the temperature equilibrium state. However, the screw shaft is under the combined effect of three heat sources, including the front bearing, the rear bearing, and the nut. Then Eq. (50) is not the exact solution of the temperature response of the screw shaft. Instead, the exact solution is the superstition of the temperature responses of three heat sources. Namely,
T ðx; tÞ ¼ T f ðx; tÞ þ T r ðx; tÞ þ T s ðx; t Þ þ T 0
ð51Þ
(5) Thermal Error Modeling and Identification Methods For the linear axis with the fixed-support or fixed-free installation, the transient thermal elongation, caused by the temperature rise, is then obtained as
RL ða0 DT ðx; t ÞÞdx ¼ 0 ða0 ½T ðx; t Þ T 0 Þdx 0 RL 0 ¼ 0 a T f ð0Þ T f ðxÞ T f ðt Þ þ T r ð0Þ T r ðxÞ T r ðt Þ þ T s ð0Þ T s ðxÞ T s ðtÞ dx
DLðt; x; T Þ ¼
RL
ð52Þ
where DT ðx; t Þ denotes the temperature rise; a0 denotes the thermal expansion coefficient; T f ð0Þ, T f ðxÞ, and T f ðt Þ denote the initial temperature, the position component, and the time component of front bearing’s temperature; T r ð0Þ, T r ðxÞ and T r ðtÞ denote the initial temperature, the position-dependent component, and the time-dependent component of rear bearing’s temperature; and T s ð0Þ, T s ðxÞ and T s ðtÞ denote the initial temperature, the position-dependent component, and the timedependent component of nut’s temperature. The Taylor expansion is conducted for both the temperature components T ðxÞ and T ðtÞ, and then they are obtained as
T ðxÞ ¼
1 X T ðkÞ ð0Þ k x k! k¼0
ð53-aÞ
T ðt Þ ¼
1 X T ðkÞ ð0Þ k t k! k¼0
ð53-aÞ
By substituting Eq. (53) into Eq. (52), it is revealed that the thermal error is approximated by the polynomials with time and position as their independent variables. Then the thermal error of the linear axis is approximated by these polynomials. It is noticed that the thermal expansion at different positions is a function with time t as its independent variable. Then it is difficult to calculate the thermal expansion of the screw shaft directly by Eq. (52) and Eq. (53). But the above analysis provides the basic idea for the thermal error fitting of the linear axis. If the high-order terms are ignored, then the thermal error of the linear axis is expressed as
Eðx; t; T Þ ¼ DLðt; x; T Þ ¼ At t 3 þ Bt t þ C t ax x3 þ bx x þ cx
ð54Þ
where At , Bt , and C t , which are related to the time, are the coefficients to be determined; and ax , bx , and cx , which are related to the position, are the coefficients to be determined. The temperature of the screw shaft is regarded as the intermediate variables to simplify the decomposition of the thermally-induced error from the measured positioning error of linear axes. The measured positioning errors of linear axes are the sum of the thermal error and the geometric error, namely,
DEðx; t; T Þ ¼ Eðx; t; T Þ þ oðxÞ
ð55Þ
where oðxÞ denotes the geometric error. The thermal error Eðx; t; T Þ is a function with time and position as its independent variables, while the geometric error oðxÞ is a function with position as its independent variable. The above operation greatly simplifies the thermal error modeling of the linear axis, and then the decomposition of the thermal error from the measured positioning error is realized. The basic idea for the modeling method of thermal error of the linear axis is then obtained.
Table 3 Parameters used in Eq. (50). Parameters
a (m/°C)
h (W/m2/ °C)
k (W/(m °C))
T 1 (°C)
T ð0Þ (°C)
L(mm)
Values
1.0e-5
70.1
60
21
60
400
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
17
Fig. 8. Temperature distribution of the linear axis.
2.3.2. Thermally-induced error modeling of spindle system The thermal error of the spindle system is closely related to temperatures of critical points because the thermal error is the manifestation of the temperature rise. The MLRA model is widely adopted because of its simple structure and highefficiency calculation. Moreover, the MLRA-based thermally-induced error model has a clear physical meaning, and then it is extremely convenient for engineering applications. The thermal error modeling method of the spindle system will be studied based on the MLRA. (1) Thermally-induced Error Modeling Based on MLRA Method The MLRA is effective to establish the correlation between the thermal error and the critical temperatures. Then the relationship between the dependent variable and the independent variable is
y ¼ b 0 þ b 1 x1 þ b 2 x2 þ þ b m xm þ e
ð56-aÞ
where b0 , b1 , b2 , bm denote the coefficients to be determined; e denotes the residue; x1 , x2 , , and xm denote the input variables; and y denotes the output variable. The least-square method is used to calculate the coefficients. The thermal error models of the spindle unit, including the thermal elongation Ez , thermal yaw cx , thermal pitch cy , the translational error gx and gy , are expressed as
0
Ez
1
0
b0E
Bc C Bb B x C B 0cx B C B B cy C ¼ B b0cy B C B Bg C B @ x A @ b0gx gy b0gy 0
b1E
b2E
b1cx
b2cx
b1cy
b2cy
b1gx
b2gx
b1gy
b1gy
10 1 1 0 1 eEz C B T1 C B bmcx C ec C CB C C B CB B xC C B T C B bmcy CB 2 C þ B ecy C C CB .. C C C B bmgx AB @ . A @ egx A egy b1gy Tm bmE
ð56-bÞ
1 0 1 b0E b1E b2E bmE eEz B b0cx b1cx b2cx bmcx C B ecx C B C B C C B C where B B b0cy b1cy b2cy bmcy C denotes the coefficient matrix; B ecy C denotes the residual matrix; and T 1 , T 2 , , and @ b0g b1g b2g bmg A @ egx A x x x x egy b0gy b1gy b1gy b1gy T m denote the critical temperatures. (2) Angular errors According to Fig. 9, the translational error of the spindle system in the Z-direction is measured by the displacement sensor S5. Moreover, the thermal yaw cy is measured by displacement sensors S1 and S3. Namely,
18
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 9. The spindle thermal inclination sketch.
tancy ¼
DY 1 DY 3 DH
ð57-aÞ
where DY 1 and DY 3 denote the translational errors of the spindle system in Y-direction; and DH ¼ 100mm. Then the thermal yaw cy is expressed as
DY 1 DY 3 DH
cx ¼ arctan
ð57-bÞ
Correspondingly, the thermal pitch cx is expressed as
tancx ¼
DX 2 DX 4 DH
ð58-aÞ
where DX 2 andDX 4 denote the translational errors of the spindle system in X-direction. Then, the thermal pitch cx is expressed as
DX 2 DX 4 DH
cx ¼ arctan
ð58-bÞ
(3) Translational errors The translational error of the spindle system in the Y-direction is mainly attributed to the angular error around Zdirection and that in Y-direction. So the translational error in Y-direction is expressed as
DX 1 ¼ gx þ R cos/ ð1 coscy Þ
ð59-aÞ
where / and R denote the rotational angle of the spindle system and the distance between two sensors on the test rod, respectively. Similarly,
DY 1 ¼ gy þ R cos/ ð1 coscx Þ
ð59-bÞ
where gy denotes the translational error in the Y-direction. So the translational errors gx and gy of the spindle system in the X- and Y-directions are then identified according to Eq. (59).
(
gx ¼ DX 1 R cos/ ð1 coscy Þ gy ¼ DY 1 R cos/ ð1 coscx Þ
ð60Þ
where DX 1 andDY 1 denote the distances in X- and Y- directions from the probe of temperature sensor S1 to the test rod, respectively.
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
19
3. Thermally-Induced error identification of linear and spindle axes 3.1. Experimental setup and measuring principles The thermal error terms, including the Eðx; t; T Þ, Eðy; t; T Þ, Eðz; t; T Þ, gx , gy , gz , cx , cy , and cz will be accurately characterized. The angular errors and the translational errors of the spindle system are measured. The angular errors of the spindle system were measured by the five-point method [43]. The displacement senor S5 was used to collect the axial thermal elongation. The displacement sensors (S1 , S3 , S2 and S4 ) were used to measure the thermal yaw cx and thermal pitch cy , as shown in Fig. 10. Then the National Instrument cRIO-9045 was used to construct the acquisition system, in which high-precision C8-2.0 displacement and PT100 sensors were used to measure the thermal errors and the temperature rise, respectively. The performance parameters of the displacement sensors are listed in Table 4. The temperature sensors are used to capture the temperatures of the main heat sources of the spindle system. The temperature sensors, which are the OMEGA PR-11-2-100-M15-(*)-E, are the PT100 sensors. The temperature sensor is installed on the principle of capturing the temperatures of the main heat sources, and then the installation positions of these temperature sensors are listed in Table 5. For the spindle system, the thermal error is closely related to the critical temperatures T 1 ~T 5 . The laser interferometer Renishaw XL80 was used to collect the positioning errors of the X-, Y-, and Z- axes. Then the thermal error will be decomposed from the measured positioning error. The positioning error of the linear axis with fixed-support installation was investigated, depicted in Fig. 11. To simulate the thermal equilibrium process, the nut with the worktable travels back and forth in the full stroke of the linear axis. Then the Renishaw laser interferometer XL80 was applied to collect the positioning error at each measurement point. The measuring principle is shown in Fig. 12. The whole stroke of the linear axis is from 0 ~ 700 mm. Then the measuring range of the Renishaw laser interferometer XL80 is set as 0~650 mm because of the installation of the laser interferometer. Here the reference point is set as X = 0 and the measuring spacing is set as 50 mm. So the total number of the measuring points is 11. For each measurement, the Renishaw laser interferometer XL80 takes about 2 s to accomplish the measurement and collection. Then the worktable with the retroreflector stagnates for 4 s to guarantee the complete the above process. The reverse overrun was set as 2 mm to avoid the measurement error caused by the backlash. The linear axis continuous operating for 9 h to reach the thermal equilibrium state and then the worktable moves back and forth for several times in the full stroke of the linear axis. The measured positioning errors are obtained as the average value of the collected data. The air pressure and humidity sensors were used to compensate for the fluctuation effect of the air pressure and humidity. The feed rate, the rotational speed, and the ambient temperature are the main factors, which have a significant effect on the thermal error of the five-axis machine tool. By changing the above three factors, then the working conditions of the fiveaxis machine tool are set according to Table 6 to simulate the actual machining process. The working conditions #1, #2, and #3 are used for the modeling, the prediction and the verification, respectively.
Y-axis
X-axis Test bar S 3 S1 S5 S4
Spindle
S2
Worktable
Z-axis Fig. 10. Diagram of displacement sensor installation.
20
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Table 4 Performance parameters of displacement sensor. Type
Manufacturer
Model
Measuring range
Bandwidth
Linearity
Peak resolution
Capacitive
Lion
C8-2.0
250 mm
10KHz
0:15% FSO
6.32 nm
Table 5 Installation of temperature sensors. Temperature sensors
Installation positions
T1 T2 T3 T4 T5
Motor stator Front bearing housing Rear bearing housing Coolant outlet Spindle Housing
5-axis machine tool
Renishaw laser interferometer XL80 PC Fig. 11. Experimental setup.
Fig. 12. Measuring principle of positioning error of ball screw feed drive system.
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
21
Fig. 13. Rotational speed distribution.
3.2. Thermally-Induced error terms identification 3.2.1. Identification of Thermally-Induced error terms of linear axis It takes 9 h to conduct the thermal characteristic experiments of the linear axis, and the measured positioning errors of X-, Y-, and Z- linear axes are then obtained, seen in Fig. 14. It is depicted that the positioning error increases with the position in the negative direction and that the increase rate declines with the position. Moreover, the positioning errors at a different time have an identical shape, which means the measurement of the positioning error is reliable. Besides, the positioning error of the linear axis increases with time because the heat conduction efficiency is higher than that of the heat convection. It is revealed that the positioning error of the linear axis shows strong nonlinearity respective to the position. The empirical model is used in previous studies to establish the correlation between the thermal error and the temperature rise, which ignores the thermal error generation mechanism of the linear axis. Then the necessity to construct the thermally-induced error correlation of linear axes from the view of the thermal error generation mechanism becomes urgent according to Section 2.3.1. For the linear axis with fixed-support installation, the thermal error is formulated as the product of the polynomial function with time as its independent variable and the polynomial function with position as its independent variable. Then the thermally-induced error correlation is used to fit the measured data of X-, Y-, and Z-axes according to Eq. (54). Here the polynomial model with degree 3 is used, and the maximum fitting goodness of R2 is 0.95, which verifies the fitting performance of the proposed model in this paper. Then the fitted positioning errors are compared with the measured data, seen in Fig. 15. It is shown that the fitted positioning errors match well with the measured data in the whole moving range and the whole thermal equilibrium process, and then the validity of the proposed method is verified. The thermal error is the main component of the positioning error of the linear axis, and then the decomposition of the thermal error is feasible according to Section 2.3.1. The geometric error oðxÞ is the difference between the measured positioning error and the fitted thermal error, and the geometric error oðxÞ at t ¼ 0 is greater than that of at other time because the thermal effect has not come into effect at t ¼ 0. To exclude the errors in the measurement and modeling, the thermal error model is modified as
22
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Table 6 Working conditions of linear axis. Working conditions
#1 #2 #3
Parameters
Functions
Feed rate of linear axis
Rotation speed of spindle
Ambient temperature
24 m/min 30 m/min 18 m/min
Fig. 13(a) Fig. 13(b) Fig. 13(c)
22:4 0:5o C 19:2 0:3o C 19:2 0:5o C
Modeling Prediction Compensation
Fig. 14. Positioning errors linear axis.
ax x3 þ bx x þ cx Eðx; t; T Þ ¼ At t 3 þ Bt t þ C t a 0 x3 þ b 0 x þ c 0
ð61Þ
where a0 , b0 , and c0 , which are related to the thermal error at x ¼ 0, denote the coefficients at x ¼ 0. Finally, the correlation between the thermally-induced error and position is established, and the thermally-induced errors of X-, Y-, and Z-axes are identified according to the above analysis. Then the identified thermal errors are then used for subsequent compensation. Then the fitting goodness of the RMSE, the fitting ability g, the absolute mean, the absolute maximum, and the absolute minimum are used to evaluate the fitting performance of these models, seen in Tables 7–9. The fitting ability g for the X-axis at 0 h, 2 h, 5 h, 7 h, and 9 h are 94.99%, 93.27%, 94.36%, 93.71%, and 93.52%, respectively. The fitting ability g for the Y-axis at 0 h, 2 h, 5 h, 7 h, and 9 h are 96.13%, 96.11%, 98.00%, 98.67%, and 98.72%, respectively. The fitting ability g for the Z-axis at 0 h, 2 h, 5 h, 7 h, and 9 h are 97.87%, 97.19%, 98.20%, 98.64%, and 98.08%, respectively. The fitting performance of the thermal error models of the feed drive systems is excellent.
23
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 15. Comparison between fitted thermal error and measured positioning errors of linear axis.
Table 7 Model evaluation of fitting performance for X-axis. Indicators
RMSE Fitting ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time(h) 0
2
5
7
9
0.0689 94.99 0.2339 0.3686 0.0129
0.1201 93.27 0.3086 0.4843 0.0150
0.2462 94.36 0.4416 0.6883 0.0185
0.3710 93.71 0.5423 0.8496 0.0259
0.5366 93.52 0.6523 1.0244 0.0328
0
2
5
7
9
0.0200 96.13 0.1273 0.2000 0.0000
0.0827 96.11 0.2636 0.4000 0.0000
0.0618 98.00 0.2364 0.3999 0.1000
0.0418 98.67 0.2000 0.3000 0.1000
0.0491 98.72 0.2182 0.3000 0.2000
Table 8 Model evaluation of fitting performance for Y-axis. Indicators
RMSE Fitting ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time(h)
24
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Table 9 Model evaluation of fitting performance for Z-axis. Indicators
RMSE Fitting ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time (h) 0
2
5
7
9
0.0282 97.87 0.1555 0.2000 0.1000
0.0964 97.19 0.2909 0.5000 0.1000
0.0618 98.20 0.2405 0.3000 0.1000
0.0564 98.64 0.2182 0.4000 0.2000
0.1382 98.08 0.3636 0.5000 0.0282
Fig. 16. Temperature rise at working condition #1.
3.2.2. Identification of Thermally-Induced error terms of spindle system The critical temperatures of the spindle system at working condition #1 are obtained, as shown in Fig. 16. The trends of temperature rise and drop are consistent with the trend of the rotational speed distribution shown in Fig. 13 (a), which demonstrates the dependence of the temperature on the rotational speed. The heat generation of bearings increases with the rotational speed according to Ref. [6]. Then the temperature is the direct reflection of the heat generation of bearings. Finally, it is concluded that the temperature rise is closely related to the rotational speed distribution of working condition #1. Moreover, the fluctuation of the temperatures means that the change in speeds causes the difficulty to achieve a heat balance. Then the thermal deformations of the high-speed spindle at working condition #1 are obtained, as shown in Fig. 17(a). The thermal deformations increase with the rotational speed of the spindle system, and there exists a strong lag effect between the thermal deformation and the rotational speed. Then the thermal yaw cx and thermal pitch cy are obtained, as shown in Fig. 17(b). Finally, the translational errors gx and gy are identified according to Eq. (60), seen in Fig. 17(c). The critical temperatures shown in Fig. 17(a) are regarded as the inputs of the MLRA model, and then the thermal errors are predicted by the critical temperatures with the aid of the MLRA method, seen in Eq. (62). The residues between the measured data and fitting errors are in the range of [-0.49 mm, 0.5 mm]. For transitional errors, the fitting goodness of R2 is greater than 0.94, which verifies the effectiveness of the modeling method.
Ez ¼ 89:6619 þ 0:5752T 1 83:8363T 2 þ 10:6994T 3 þ 14:1400T 4 12:6225T 5
cx ¼ 31:1984 0:5500T 1 12:6472T 2 þ 7:8383T 3 þ 9:6489T 4 3:0353T 5 cy ¼ 78:1968 þ 2:6834T 1 þ 22:0344T 2 21:7061T 3 18:7970T 4 þ 12:4126T 5 gx ¼ 177:2860 0:5828T 1 11:2493T 2 þ 11:5722T 3 þ 4:3105T 4 þ 3:8117T 5 gy ¼ 72:7147 0:2376T 1 4:5866T 2 þ 4:7182T 3 þ 1:7575T 4 þ 1:5541T 5
ð62Þ
Then Table 10 lists the fitting goodness for the evaluation of the predictive performance. The fitting ability g of the models for Ez , cx , cy , gx , and gy are 96.28%, 96.10%, 96.44%, 94.47%, and 89.10%, respectively. The fitting performance of the thermal error models of the spindle system is excellent. 4. Modal validation 4.1. Thermally-Induced error prediction 4.1.1. Prediction of Thermally-Induced error terms of linear axis To verify the validity of the thermal error modeling method, the static thermal error models established in Section 3.2.1 are then used to predict the thermal errors of the linear axis at working condition #2, and the predicted thermal errors are
25
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 17. Thermal error of the spindle system at working condition #1.
Table 10 Model evaluation of fitting performance. Indicators
E (lm)
cx (00 )
cy (00 )
gx (lm)
gy (lm)
RMSE Fitting ability g (%) Mean of absolute Absolute maximum Absolute minimum
0.5815 96.28 0.7148 1.0417 0
0.5783 96.10 0.7114 1.0523 0
0.5802 96.44 0.7137 1.0382 0
2.9336 94.47 1.3662 4.5612 0.0204
4.1275 89.10 1.6662 4.3602 0.0274
compared with the measured data, as shown in Fig. 18. The residues between the predicted errors and the measured data are in the range of [-0.6 mm, 0.6 mm]. The predicted errors match well with the measured data in the whole measuring range, which validates the effectiveness and robustness of the thermal error modeling method. Then the prediction goodness of the RMSE, the fitting ability g, the absolute mean, the absolute maximum, and the absolute minimum are used to evaluate the fitting performance of these models, seen in Tables 11–13. The predictive ability g of these models has no significant decrease, which validates that the models are robust and generalizable. The excellent predictive performance and strong robustness come from the mechanism error modeling of the feed drive system.
4.1.2. Prediction of thermally-induced error terms of spindle system The thermal characteristic of the spindle system at working condition #2, including the temperature rise and the thermal deformations, was measured, seen in Fig. 19. It is depicted that both the temperature and thermal deformations show a strong nonlinear relationship with the running time. Moreover, the dependence of the temperature and thermal deformations on the rotational speed distribution is still significant. The maximum axial thermal elongation reaches about 50 mm, and the maximum radial deformation is about 40 mm.
26
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 18. Comparison between fitted thermal error and measured positioning errors of linear axis at working condition #2.
Table 11 Model evaluation of prediction performance for X-axis. Indicators
RMSE Predictive ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time(h) 0
2
5
7
9
0.4144 88.57 0.5649 1.0825 0.1760
0.5703 90.24 0.6287 1.2870 0.1014
0.5502 93.41 0.6250 1.2684 0.1000
0.3571 95.52 0.5370 1.1765 0.1512
0.7028 94.83 0.6592 2.0923 0.0127
The thermal errors of the spindle system are then predicted with the critical temperatures in Fig. 19(a) as the inputs. The predicted thermal errors and the measured data are then compared, seen in Fig. 20. It is revealed that the thermal error models established in Section 3.2.2 can predict the thermal errors of the spindle system accurately, which validates the effectiveness and robustness of the proposed method. Then Table 14 lists the predictive goodness for the evaluation of the predictive performance. The predictive ability g of the models for E, cx , cy , gx , and gy are 98.42%, 92.10%, 95.97%, 95.90%, and 96.27% respectively. The predictive performance of the thermal error models is great. Comparing with the fitting ability g, then the robustness and generalizable ability are validated.
27
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538 Table 12 Model evaluation of prediction performance for Y-axis. Indicators
RMSE Predictive ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time(h) 0
2
5
7
9
0.0355 94.41 0.1495 0.3366 0.0074
0.0269 96.74 0.1230 0.3837 0.0162
0.0236 97.86 0.1198 0.3215 0.0147
0.0367 98.28 0.1437 0.3432 0.0156
0.0291 98.83 0.1438 0.2722 0.0077
0
2
5
7
9
0.1053 92.01 0.2666 0.5190 0.0136
0.0509 96.22 0.1780 0.4500 0.0247
0.5094 91.02 0.6084 1.1661 0.1310
0.0363 98.66 0.1349 0.3473 0.0033
0.0634 98.49 0.2066 0.4478 0.0231
Table 13 Model evaluation of prediction performance for Z-axis. Indicators
RMSE Predictive ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time(h)
Fig. 19. Thermal characteristics of spindle system at working condition #2.
4.2. Thermally-induced error compensation To improve the machining accuracy of five-axis machine tools, the thermal error compensation experiments were conducted on the CNC machine tools. The thermal error compensation principle of the spindle system is shown in Fig. 21. The thermal error compensation system includes the thermal characteristic acquisition and error compensation modules. The PT100 temperature sensors are used to collect the temperatures. The running time and actual position are collected from the CNC system. Then these signals are filtered, amplified, and A/D converted to get the actual physic information of machine tools. Then the thermally-induced error terms are identified. For the spindle system, the thermally-induced error is identified as a function with critical temperatures as its independent variables. For the linear axes, the thermal errors are identified as the product of the polynomial function with time as its independent variable and the polynomial function with position as its independent variable. Finally, the thermal error compensation components are calculated in each direction according to Eq. (23). The compensation components of the axis are sent to the CNC controller by the compensation module. The reverse of the error is superimposed to the output commands of the servo controller to offset the position deviation by comparing the actual and command positions of the worktable. Finally, the real-time thermal error compensation is realized. The five-axis machine tool is warmed up at working condition #3. Then thermal error compensation is carried out for linear axes. Moreover, the positioning errors of the linear axes are measured separately. The compensation results are then compared for X-, Y-, and Z- axes, as shown in Fig. 22. It is shown that the present error compensation can effectively reduce
28
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 20. Thermal error of spindle system at working condition #2.
Table 14 Model evaluation of prediction performance. Indicators
E (lm)
cx (00 )
cy (00 )
gx (lm)
gy (lm)
RMSE Predictive ability g (%) Mean of absolute Mean of absolute maximum Mean of absolute minimum
0.3234 98.42 0.4836 0.9937 0.0289
1.0561 92.10 0.9412 1.9351 0.1216
0.3205 95.97 0.4966 0.9441 0.0046
1.2083 95.90 0.9403 1.9589 0.0086
0.5946 96.27 0.6602 1.4273 0.0000
Fig. 21. Schematic diagram of compensation for thermal error of spindle system.
the positioning error, and the traditional pitch error compensation can also reduce the positioning error to some extent, which verifies the validity of the proposed method. The positioning error with the present error compensation is much smaller than that with the traditional pitch error compensation, which indicates that the necessity for the consideration of dynamic changes in the thermal error. For the present error compensation, the positioning errors of X-, Y-, and Z-axes are in the range of [4.8 lm, 3.5 lm], [5.5 lm, 3.0 lm], and [4.8 lm, 1.9 lm], respectively. For the traditional pitch error
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
29
Fig. 22. Comparison of error compensation at working condition #3.
compensation, the positioning errors of X-, Y-, and Z-axes are in the range of [9.8 lm, 9.7 lm], [10.9 lm, 8.1 lm], and [11.8 lm, 6.2 lm], respectively. The positioning errors of X-, Y-, and Z-axes are in the range of [4.8 lm, 3.5 lm], [5.5 lm, 3.0 lm], and [4.8 lm, 1.9 lm], respectively. For the linear axes with thermal error compensation, the positioning error is reduced significantly, which verifies the effectiveness of the thermal error compensation method and validates the statement that the thermal error is the main component of the positioning error further. It is concluded that the positioning errors of the linear axes are reduced obviously and that the thermally-induced error compensation method is effective and robust. Then the thermally-induced error compensation is implemented for the spindle system at working condition #3. The compensation results are then obtained, as shown in Fig. 23. It is shown that the thermally-induced errors of the spindle system are decreased significantly. The thermal elongation, the thermal yaw cx , and thermal pitch cy are in the range of [3.4 lm, 5.0 lm], [3.000 , 2.500 ], and [3.000 , 3.000 ], respectively. It is noticed that the thermal errors start to decreases from about 275 min, which is slightly lagging behind the decrease in the rotational speed shown in Fig. 13(c). 4.3. Actual Machining The validity of the data-driven thermal error compensation method is fully validated when actual machining is performed. So a large number of test samples were machined by the five-axis machine tool with and without thermallyinduced error compensation at working condition #3. Before the actual machining of test specimens, the motion simulation is conducted for the five-axis machine tool to simulate the actual machining process. The SS304 is chosen as the test material because of its wide use, as shown in Fig. 24. The workpiece is fixed on the B-axis and then the linkage of the X-, Y-, Z-, C-, B-, and S- axes realize the precision machining of workpiece. The feed rates of X-, Y-, and Z- axes are 1.2 m/min. The rotational speeds of the C-axis and B-axis are 0. The rotational speed of the spindle system is 5000 r/min. The cutting depth is 50 mm. It takes about one day to machine 32 test specimens. When the machining of all test specimens is finished, these specimens are
30
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 23. Comparison of the compensation results at working condition #3.
D2 D1
Fig. 24. Test samples.
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
31
Fig. 25. Machining error of the test samples.
placed in the laboratory with a constant temperature of 20 °C for 8 h, and then the roundness error of the test specimens is measured by the coordinate measuring machine, as shown in Fig. 25. The sizes of D1 and D2 are the most important parameters for the test specimens because of the assembly requirement. The basic dimensions of D1 and D2 are 70 mm and 95 mm, respectively. Then the roundness error of the workpiece is measured with the coordinate measuring machine. It is seen that the roundness error distribution is much closer to the zero line and that the roundness error is reduced by more than 85% and 37% with the present error compensation compared with that without compensation and that without traditional error compensation, respectively. Then the effectiveness of the thermally-induced compensation method is totally verified. Moreover, the trend of the roundness error without compensation is similar to the preload relaxation process caused by the temperature rise of the linear axis. 5. Results and discussion 5.1. Effect of heat source on thermal error Then the temperature of the moving nut and the bearings are obtained as
T n ¼ cnQm2 n þ T 0
ð63Þ
T f 1 ¼ T f 2 ¼ c Qm1 þ T 0 b
b
where T f 1 and T f 2 denote the temperatures of the bearing #1 and bearing #2, respectively; T n denotes the temperature of the moving nut; cb and mb denote the specific heat capacity and mass of the bearing, respectively; and cn and mn denote the specific heat capacity and mass of the nut, respectively. According to Eq. (24), Eq. (35), and Eq. (63), the temperatures of the moving nut and the bearings are obtained. Two ball screws, whose models are W4009SA-6D-C5Z12, are used in the dual-drive servo feed drive system. The material of the screw shaft is ASTM E52100. The density q, the specific heat capacity c, and the thermal expansion coefficient k are obtained. The diameter d0 of the screw shaft is 85 mm. Then the forced convective coefficient h is calculated according to Eq. (38), as listed in Table 2. Then the temperature of the screw shaft is calculated, seen in Fig. 26. The measured data are
Fig. 26. Temperature response of screw shaft.
32
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 27. Positioning error of feed drive system.
compared with the fitting results. It is shown that the calculated temperature agree well with the measured data and that the temperature response of the screw shaft increases with the running time t, which is consistent with the results obtained by Ref. [44]. The temperature rise T ðx; tÞ of the screw shaft is obtained. Then the temperature rise T ðx; tÞ is substituted into Eq. (52) to calculate the thermal error of the feed drive system, seen in Fig. 27. It is seen that the calculated thermal error agrees well with the scatter points, which validates the effectiveness of the thermal error model. It is concluded that the effect of the heat source on the thermal error is realized by changing the temperatures of the heat sources. Then the heat entering linear axes is changed. By substituting the change in the heat generation into the thermal field and thermal error models of linear axes, the temperature and thermal error are predicted accurately. 5.2. Effect of pre-stretch of lead screw on thermal error When the ball screw feed system is working, a large amount of heat is generated, and the heat causes the thermal expansion of the screw shaft, and hence the lead is increased, leading to the decrease in the positioning accuracy. To compensate the thermal expansion, the pre-stretch of the ball screw is used as the compensation method. At the time of manufacturing, the length of the threaded portion of the ball screw is less than the nominal length (lead by the number of turns) a pre-stretch. The pre-stretch capacity is slightly larger than the thermal expansion. In the assembly process, the lead screw is stretched by a certain tensile structure so that the threaded portion of the lead screw reaches the nominal length. In the service process, the thermal expansion counteracts a part of the pre-stretch capacity, and then the tensile stress of the lead screw decreases, but the length of the screw shaft does not change. The above method ensures that the pitch accuracy is not affected by the thermal expansion. Under no-load conditions, the preload F t applied by the certain tensile structure is equal to the axial load F a . So F a ¼ F t . The greater the preload F t , the more is the heat generation of the screw shaft. Moreover, the greater the pre-stretch capacity, the greater is the preload F t . So the heat generation of the screw shaft increases with the pre-stretch capacity. The change in the pre-stretch capacity leads to the change in the heat generation of the screw shaft and the pre-stretch capacity and the heat generation of the screw shaft are positively correlated. The ball screw is fixed at one end and an axial force is applied to the other end. When the temperature rises, the lead screw generates the compressive stress. The thermal stress distributes in the longitudinal direction, and the induced elongation is expressed as
DL2 ¼ Lr=E
ð64Þ
where r denotes the stress and r ¼ P=A, wherein A and P denote the cross-section area and the axial force, respectively. So the actual elongation is expressed as
Z
L
DL1 ¼ DLðx; tÞ þ DL2 ¼ DLðx; t Þ ¼
ða½T ðx; tÞ T 0 Þdx Lr=E
ð65Þ
0
The pre-stretch capacity is expressed as
DL3 ¼ a D T L
ð66Þ
where D T denotes the average temperature rise. Under no-load conditions, the preload F t applied by the preload nut is equal to the axial load F a . So F a ¼ F t . The positioning error of the feed drive system is expressed as Eq. (67). To ensure that the ball screw feed drive system works properly, it
33
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538 Table 15 Temperature rises and thermal errors. Temperature rise(°C)
Thermal errors(mm)
Pre-stretch capacities (°)
Moving nut
Front bearing housing
Rear rearing housing
Maximum
0 120 150 180
3.92 3.88 4.04 3.80
4.48 4.64 5.11 5.86
3.79 3.91 4.07 5.00
43.7 41.1 32.3 30.2
is necessary to ensure that the difference between the actual elongation of the screw shaft and the pre-stretch capacity is not less than zero. Namely,
DL ¼ DL1 DL3 P 0
ð67Þ
According to Eq. (67), the increase in the pre-stretch capacity leads to the decrease in the positioning error of the linear axis with the fixed-support installation. The pre-stretch capacities are 0°, 120°, 150°, and 180° and then the thermal characteristic tests were carried out. Table 15 shows the temperature rises of the critical points and the maximum thermal errors under various pre-stretch capacities. It is seen that the larger the pre-stretch capacity, the smaller is the maximum thermal error of the feed drive system, and the better is the precision stability of the machine tool. According to Eq. (67), the positioning error of the feed drive system will be reduced with the increase in the pre-stretch capacity, which is consistent with that shown in Table 15. However, it should be noted that the pre-stretch capacities should not be as large as possible, and the excessive pre-stretch capacity will cause the excessive temperature rise of bearings. If the temperature rise is too high, the wear will intensify and even seizure happens to the bearing. Moreover, it is seen that the larger the pre-stretch capacity, the larger the temperature rise of the bearing housing because the increase in the pre-stretch capacity leads to an increase in the bearing load, and the increase in the friction heat causes the increase in the temperature of the bearing housing. In different pre-stretch states, the temperature rise of the nut does not change much because the pre-stretch capacity of the lead screw is still relatively small, and then the effect of the pre-stretch capacity on the fit between the nut and the screw is small. It is revealed that the effect of the pre-stretch of the lead screw on the thermal error of linear axes is realized by changing the deformation of the screw shaft. The pre-stretch can compensate for part of the thermal elongation of the screw shaft. The effect of the pre-stretch of the lead screw on the temperatures and thermal error of linear axes is significant. Moreover, the greater the pre-stretch of the lead screw, the higher is the temperature of linear axes, and the smaller the thermal error of linear axes. To balance the temperature and thermal error, a reasonable pre-stretch capacity should be chosen. 6. Conclusions The data-driven thermally-induced error compensation method was proposed for the five-axis machine tools based on the homogeneous transformation method. The compensation component was obtained by analyzing the error transmission chain of five-axis machine tools and was expressed as the matrix of the error terms of the linear axes and spindle system. From the view of error generation mechanism, the thermally-induced error of the linear axis was formulated as the product of the polynomial function with time as its independent variable and the polynomial function with position as its independent variable. Then the thermal and geometric errors were decomposed and identified from the measured positioning error. The thermal error of the spindle system was expressed as the MLRA model of the critical temperatures. Finally, the error compensation was carried to verify the effectiveness of the data-driven thermally-induced error compensation method. The following conclusions are obtained from the previous analysis. (1) The systematic data-driven thermally-induced error compensation method of the five-axis machine tools is proposed based on the homogeneous transformation method. Then the relationship between the compensation components and the thermal error components in the linear axes and spindle system is revealed by the differential movement of the compensation axis. The thermally-induced error model is robust and effective at different thermal effects, and the data-driven thermally-induced error compensation method could reduce by the machining error more than 85% and 37% with the present error compensation compared with that without compensation and that without traditional error compensation, respectively. (2) The relationship between the thermal expansion, the thermal error, and the resultant positioning error of the linear axis is explored. The temperature of the screw shaft is regarded as the intermediate variables to simplify the decomposition of the thermally-induced error from the measured positioning error of linear axes. For the linear axis with fixed-support installation, the thermal error is formulated as the product of the polynomial function with time as its independent variable and the polynomial function with position as its independent variable. With the aid of the
34
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
polynomial function, the thermal error modeling of the linear axis becomes ever simple and easy. The thermal error is the main component of the positioning error of the linear axis. The decomposition of the thermal error from the positioning error is realized by the polynomial fitting. (3) The direct and accurate calculation of the thermal error of the linear axis is feasible from the temperature rise. The temperature response of the screw shaft is the superposition of three responses of equivalent temperatures. The thermal error is the integral of the temperature rise on the length of the lead screw. It is revealed that the MLRA method is effective to establish the correlation between the thermally-induced error and critical temperatures of the spindle system. The dependence of the thermal error on the feed rate and rotational speed is strong. Both the temperature rise and thermally-induced error of the spindle system increase with the rotational speed and feed rate. The lag effect of the thermally-induced error is obvious respect to the fast change in the feed rate and rotational speed. (4) The effect of the heat source on the thermal error is realized by changing the temperatures of the heat sources. Then the heat entering linear axes is changed. By substituting the change in the heat generation into the thermal field and thermal error models of linear axes, the temperature and thermal error are predicted accurately. The effect of the prestretch of the lead screw on the thermal error of linear axes is realized by changing the deformation of the screw shaft. The pre-stretch can compensate for part of the thermal elongation of the screw shaft. The effect of the pre-stretch of the lead screw on the temperatures and thermal error of linear axes is significant. Moreover, the greater the prestretch of the lead screw, the higher is the temperature of linear axes, and the smaller the thermal error of linear axes. To balance the temperature and thermal error, a reasonable pre-stretch capacity should be chosen.
Acknowledgments This research was supported by the National Natural Science Foundation of China (51905057), the Natural Science Foundation Project of Chongqing, Chongqing Science and Technology Commission (cstc2019jcyj-msxmX0050), and the Fundamental Research Funds for the Central Universities (2018CDXYJX0019). Declaration of Competing Interest We declare that we have no financial and personal relationships with other people or organizations that can inappropriately influence our work; there is no professional or other personal interest of any nature or kind in any product or company that could be construed as influencing the position presented in, or the review of, the manuscript. References [1] J. Mayr, J. Jedrzejewski, E. Uhlmann, et al, Thermal issues in machine tools, CIRP Ann.-Manuf. Technol. 61 (2) (2012) 771–791. [2] Y. Altintas, M.R. Khoshdarregi, Contour error control of CNC machine tools with vibration avoidance, CIRP Ann.-Manuf. Technol. 61 (1) (2012) 335–338. [3] C. Ma, X. Mei, J. Yang, et al, Thermal characteristics analysis and experimental study on the high-speed spindle system, Int. J. Adv. Manuf. Technol. 79 (1–4) (2015) 469–489. [4] C. Ma, L. Zhao, X. Mei, et al, Experimental and simulation study on the thermal characteristics of the high-speed spindle system, ARCHIVE Proc. Inst. Mech. Eng. Part C J. Mech. Eng. Sci. 1989-1996 203–210 (6) (2016) 6. [5] C. Ma, J. Yang, L. Zhao, et al, Simulation and experimental study on the thermally induced deformations of high-speed spindle system, Appl. Therm. Eng. 86 (2015) 251–268. [6] J. Liu, C. Ma, S. Wang, et al, Thermal-structure interaction characteristics of a high-speed spindle-bearing system, Int. J. Mach. Tools Manuf. (2019), https://doi.org/10.1016/j.ijmachtools.2018.10.004. [7] J. Liu, C. Ma, S. Wang, et al, Thermal contact resistance between bearing inner ring and shaft journal, Int. J. Therm. Sci. 138 (2019) 521–535. [8] J. Liu, C. Ma, S. Wang, Thermal contact conductance between rollers and bearing rings, Int. J. Therm. Sci. (2020,), https://doi.org/10.1016/j. ijthermalsci.2019.106140. [9] C. Ma, L. Zhao, H. Shi, et al, A geometrical-mechanical-thermal predictive model for thermal contact conductance in vacuum environment, Proc. Inst. Mech. Eng., Part B: J. Eng. Manuf. 230 (8) (2016) 1451–1464. [10] A.L. Sweet, J.F. Tu, Tolerance design for the fit between bore and shaft for precision assemblies with significant error-scaling problems, Int. J. Prod. Res. 45 (22) (2007) 5223–5241. [11] C. Xia, J. Fu, J. Lai, et al, Conjugate heat transfer in fractal tree-like channels network heat sink for high-speed motorized spindle cooling, Appl. Therm. Eng. 90 (2015) 1032–1042. [12] J. Liu, C. Ma, S. Wang, et al, Contact stiffness of spindle-tool holder based on fractal theory and multi-scale contact mechanics model, Mech. Syst. Sig. Process. 119 (2019) 363–379. [13] J. Liu, C. Ma, S. Wang, Precision loss modeling method of ball screw pair, Mech. Syst. Sig. Process. (2020,), https://doi.org/10.1016/j.ymssp.2019.106397. [14] E. Creighton, A. Honegger, A. Tulsian, et al, Analysis of thermal errors in a high-speed micro-milling spindle, Int. J. Mach. Tools Manuf. 50 (4) (2010) 386–393. [15] H. Zhao, J. Yang, J. Shen, Simulation of thermal behavior of a CNC machine tool spindle, Int. J. Mach. Tools Manuf. 47 (6) (2007) 1003–1010. [16] J.J. Kim, Y.H. Jeong, D.W. Cho, Thermal behavior of a machine tool equipped with linear motors, Int. J. Mach. Tools Manuf. 44 (7) (2004) 749–758. [17] B. Bossmanns, J.F. Tu, A power flow model for high speed motorized spindles—heat generation characterization, J. Manuf. Sci. Eng. 123 (3) (2001) 494– 505. [18] B. Bossmanns, J.F. Tu, A thermal model for high speed motorized spindles, Int. J. Mach. Tools Manuf. 39 (9) (1999) 1345–1366. [19] C.H. Chien, J.Y. Jang, 3-D numerical and experimental analysis of a built-in motorized high-speed spindle with helical water cooling channel, Appl. Therm. Eng. 28 (17) (2008) 2327–2336. [20] J. Liu, C. Ma, S. Wang, et al, Thermal boundary condition optimization of ball screw feed drive system based on response surface analysis, Mech. Syst. Sig. Process. 121 (2019) 471–496. [21] Y. Kang, C.W. Chang, Y. Huang, et al, Modification of a neural network utilizing hybrid filters for the compensation of thermal deformation in machine tools, Int. J. Mach. Tools Manuf. 47 (2) (2007) 376–387.
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
35
[22] J.P. Kruth, P. Vanherck, V.D.B. Christophe, Compensation of static and transient thermal errors on CMMs, Ann. CIRP 50 (1) (2001) 377–380. [23] H. Yang, J. Ni, Dynamic neural network modeling for nonlinear, nonstationary machine tool thermally induced error, Int. J. Mach. Tools Manuf. 45 (4) (2005) 455–465. [24] R. Ramesh, M.A. Mannan, A.N. Poo, Error compensation in machine tools — a review: Part II: thermal errors, Int. J. Mach. Tools Manuf. 40 (9) (2000) 1257–1284. [25] S.H. Yang, K.H. Kim, K.P. Yong, Measurement of spindle thermal errors in machine tool using hemispherical ball bar test, Int. J. Mach. Tools Manuf. 44 (2–3) (2004) 333–340. [26] W.T. Lei, Y.Y. Hsu, Accuracy enhancement of five-axis CNC machines through real-time error compensation, Int. J. Mach. Tools Manuf. 43 (9) (2003) 871–877. [27] Y. Zhang, J. Yang, H. Jiang, Machine tool thermal error modeling and prediction by grey neural network, Int. J. Adv. Manuf. Technol. 59 (9–12) (2012) 1065–1072. [28] J. Yang, J. Yuan, J. Ni, Thermal error mode analysis and robust modeling for error compensation on a CNC turning center, Int. J. Mach. Tools Manuf. 39 (9) (1999) 1367–1381. [29] Y. Hong, J. Ni, Adaptive model estimation of machine-tool thermal errors based on recursive dynamic modeling strategy, Int. J. Mach. Tools Manuf. 45 (1) (2005) 1–11. [30] C. Ma, L. Zhao, X. Mei, et al, Thermal error compensation based on genetic algorithm and artificial neural network of the shaft in the high-speed spindle system, Proc. Inst. Mech. Eng., Part B: J. Eng. Manuf. 231 (5) (2017) 753–767. [31] C. Ma, L. Zhao, X. Mei, et al, Thermal error compensation of high-speed spindle system based on a modified BP neural network, Int. J. Adv. Manuf. Technol. 89 (9–12) (2017) 3071–3085. [32] Z.Z. Xu, C. Chang, L.J. Liang, et al, Study on a novel thermal error compensation system for high-precision ball screw feed drive (2nd report: Experimental verification), Int. J. Precis. Eng. Manuf. 16 (10) (2015) 2139–2145. [33] W.S. Yun, S.K. Kim, W.C. Dong, Thermal error analysis for a CNC lathe feed drive system, Int. J. Mach. Tools Manuf. 39 (7) (1999) 1087–1101. [34] H. Shi, B. He, Y. Yue, et al, Cooling effect and temperature control of oil cooling system for ball screw feed drive system of precision machine tool, Appl Therm Eng. 161 (10) (2019) 114150. [35] Z.Z. Xu, X.J. Liu, S.K. Lyu, Study on positioning accuracy of nut/shaft air cooling ball screw for high-precision feed drive, Int. J. Precis. Eng. Manuf. 15 (1) (2014) 111–116. [36] H. Pahk, S.W. Lee, Thermal error measurement and real time compensation system for the CNC machine tools incorporating the spindle thermal error and the feed axis thermal error, Int. J. Adv. Manuf. Technol. 20 (7) (2002) 487–494. [37] A.C. Okafor, ERTEKIN, et al, Derivation of machine tool error models and error compensation procedure for three axes vertical machining center using rigid body kinematics, Int. J. Mach. Tools Manuf. 40 (8) (2000) 1199–1213. [38] Z. Li, J. Yang, K. Fan, et al, Integrated geometric and thermal error modeling and compensation for vertical machining centers, Int. J. Adv. Manuf. Technol. 76 (5–8) (2015) 1139–1150. [39] Y. Wang, G. Zhang, K.S. Moon, et al, Compensation for the thermal error of a multi-axis machining center, J. Mater. Process. Technol. 75 (1–3) (1998) 45–53. [40] J.H. Jung, J.P. Choi, S.J. Lee, Machining accuracy enhancement by compensating for volumetric errors of a machine tool and on-machine measurement, J. Mater. Process. Technol. 174 (1–3) (2006) 56–66. [41] Y. Lu, A new approach to thermally induced volumetric error compensation, Int. J. Adv. Manuf. Technol. 62 (9–12) (2012) 1071–1085. [42] Incropera, Frank P. Fundamentals of Heat and Mass Transfer. Wiley, 1985:139–162. [43] ISO 230-3, Test Code for Machine Tools Part 3: Determination of Thermal Effects, ISO copyright office, Switzerland, 2007.. [44] H. Shi, C. Ma, J. Yang, et al, Investigation into effect of thermal expansion on thermally induced error of ball screw feed drive system of precision machine tools, Int. J. Mach. Tools Manuf. 97 (2015) 60–71.
Contents lists available at ScienceDirect
Mechanical Systems and Signal Processing journal homepage: www.elsevier.com/locate/ymssp
Data-driven thermally-induced error compensation method of high-speed and precision five-axis machine tools Jialan Liu a,b, Chi Ma a,b,⇑, Shilong Wang a,b a b
College of Mechanical Engineering, Chongqing University, Chongqing 400044, China State Key Laboratory of Mechanical Transmission, Chongqing University, Chongqing 400044, China
a r t i c l e
i n f o
Article history: Received 7 September 2019 Received in revised form 8 November 2019 Accepted 22 November 2019
Keywords: Five-axis machine tool Thermally-induced error Homogeneous transformation Error compensation
a b s t r a c t The data-driven thermal error compensation method of high-speed and precision five-axis machine tools was proposed based on the homogeneous transformation. The compensation component was obtained by analyzing the error transmission chain of machine tools and was expressed as the transmission matrix consisting of thermal error terms of linear axes and spindle system according to the differential movement of the compensation axis. From the view of thermal error generation mechanism, the thermal error of the linear axis was formulated as the product of the polynomial function with time as its independent variable and the polynomial function with position as its independent variable. Then the thermal errors of the linear axes were decomposed and identified from the measured positioning error. The thermal error of the spindle system was identified by the thermal characteristic experiments. The compensation component in each direction is calculated by substituting the identified thermal error terms of the linear axes and spindle system into the data-driven thermal error model. Finally, to demonstrate the effectiveness of the proposed method, the thermal error at a new working condition was predicted, and then the error compensation and actual machining were carried out. The results show that the machining error is reduced by more than 85% and 37% with the present error compensation compared with that without compensation and that without traditional error compensation, respectively. This research sheds new light on both the generation mechanism and the compensation method of the thermally-induced error of five-axis machine tools. Ó 2019 Elsevier Ltd. All rights reserved.
1. Introduction The demand of modern manufacturing to five-axis CNC machine tools is high-precision and high-speed machining, resulting from the pursuit of the machining quality and efficiency. High-speed and precision five-axis CNC machines tools are used to machine complex and precision parts in the automotive, aerospace, and military industry because it enables complex spatial movement between the tool and the workpiece. In particular, the use of high-speed and precision five-axis machine tools improves the machining accuracy and efficiency greatly. However, the high-speed and precision machine tools are under the combined effect of the geometric, dynamic, and thermal errors. Moreover, the thermally-induced error of the precision machine tools, caused by the uneven temperature distribution, accounts for up to 75% of the total volumetric error of machined workpieces [1]. So the precision of parts will be increased significantly if the thermally-induced error of ⇑ Corresponding author. E-mail addresses: [email protected] (J. Liu), [email protected] (C. Ma), [email protected] (S. Wang). https://doi.org/10.1016/j.ymssp.2019.106538 0888-3270/Ó 2019 Elsevier Ltd. All rights reserved.
2
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Nomenclature O1 O2 O3 O4 x y z b c s t B C S X Y Z
the intersection of B- and C-axes origin of the coordinate system b origin of the coordinate system c origin of the coordinate system s Local coordinate system fixed to X-axis Local coordinate system fixed to Y-axis Local coordinate system fixed to Z-axis Local coordinate system fixed to B-axis Local coordinate system fixed to C-axis Local coordinate system fixed to spindle system Coordinate system of tool nose B-axis C-axis S-axis X-axis Y-axis Z-axis a rotational angle of C-axis b rotational angle of B-axis kx angular errors in the compensation movement process ky angular errors in the compensation movement process kz angular errors in the compensation movement process rx translational errors in the compensation movement process ry translational errors in the compensation movement process rz translational errors in the compensation movement process dw change in coordinate system caused by differential transformation t T M2 friction torque between nut and screw shaft rotational speed of the screw shaft n2 Ph lead of screw shaft; V feed rate b0 helix angle of raceway Me friction torque caused by the external load Mg friction torque caused by geometric slip of balls mb and ma coefficients related to eccentricity of contact deformation ellipse # parameter related to material properties P q sum principle curvatures M gy gyroscopic moment F cq centrifugal force, F iq and F eq shear loads Q iq and Q eq normal loads A1 and A2 x and z coordinates of inner and outer groove curvature centers M1 friction torque of the bearings friction torque caused by the external loads Ml a the contact angle of the bearing. Max the greater value between 0:9F a =tana 0:1F r and F r Fa and Fr bearing’s axial and radial loads vo lubricant kinematic viscosity at operating temperature N rotational speed of bearing fo constant related to the bearing model and lubrication type Re Reynolds number Pr Prandtl number w angular velocity v viscosity of air diameter of the screw shaft d0 n1 rotational speed of screw shaft x and t the position and the time L length of screw shaft TL temperature of screw shaft C 1 and C 2 Coefficients related to temperatures
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
T ð xÞ position-dependent of temperature T ðt Þ time-dependent components of temperature T r ð xÞ position-component of rear bearing’s temperature T s ð0Þ initial temperature of nut’s temperature time-component of nut’s temperature T s ðt Þ oðxÞ geometric error Eðx; t; T Þ thermal error b0 , b1 , b2 , bm coefficients to be determined y output variable DY 1 and DY 3 translational errors in Y-direction cx thermal pitch / rotational angle of spindle system a0 , b0 , and c0 coefficients at x ¼ 0 gx and gy translational errors T f 1 and T f 2 temperatures of bearing #1 and bearing #2, specific heat capacity and mass of bearing cb and mb Ft preload r stress DL2 induced elongation DL3 pre-stretch capacity a0 thermal expansion coefficient y moving distance of Y-axis d1 offset of O2 in the x-direction of the reference system O1 d2 offset of O3 in the y-direction of the reference system O1 eYX squareness error between Y- and X-axes; Eðy; t; T Þ thermally-induced error of Y-axis dxy ; dyy ; dzy ; exy ; eyy and ezy six position-dependent geometric errors of Y-axis x moving distance of X-axis dxx ; dyx ; dzx ; exx ; eyx and ezx six position-dependent geometric errors of X-axis Eðx; t; T Þ thermal error of X-axis z moving distance of Z-axis dxz ; dyz ; dzz ; exz ; eyz and ezz six position dependent geometric errors of Z-axis eZX squareness error between Z- and X- axes eZY squareness error between Z- and Y- axes Eðz; t; T Þ thermal error of Z-axis L length of cutter d3 offset of O4 in z- direction of the reference system O1 E identity matrix cx , cy , and cz angular errors in x-, y-, and z- directions gx , gy , and gz translational errors in x-, y-, and z- directions ð X; cx Þ transformation matrixes around three x-axis R transformation matrixes around three y-axis R Y; cy R Z; cz transformation matrixes around three z-axis w homogeneous transformation matrix of the compensation movement of coordinate systems t to w t Tc Dx differential compensation movement of X-axis Dy differential compensation movement of Y-axis Dz differential compensation movement of Z-axis Da differential compensation movement of C-axis Db differential compensation movement of B-axis B the inner and outer groove curvature constant; a0 the initial contact angle D the diameter of balls h the angular displacement of the inner ring wq position angle of qth ball Ri the radius of the groove curvature center contact deformations of balls and grooves da and dr K i and K o contact stiffness aeq and aiq contact angles Fa axial load radial load Fr Mg moment load
3
4
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
X 1 and X 2 initial values of x and z coordinates of balls n1 rotational speed of the screw shaft Mv friction torque caused by lubricant viscous friction f1 coefficient related to the bearing structure equivalent load applied on bearings Fb dm pitch diameter of bearings L feature size k thermal conductivity of fluid Nu Nusselt number T0 temperature of ambient; h convective coefficient q density and the respectively; c specific heat capacity d0 denote diameter of screw shaft T0 reference temperature k thermal conductivity T ð0Þ temperature of heat source DT ðx; t Þ temperature rise initial temperature of front bearing’s temperature T f ð0Þ T f ð xÞ position-component of front bearing’s temperature T f ðt Þ time-component of front bearing’s temperature T r ð0Þ initial temperature of rear bearing’s temperature time-component of rear bearing’s temperature T r ðt Þ T s ðxÞ position-component of nut’s temperature At , Bt , and C t coefficients related to time ax , bx , and cx coefficients related position e the residue x1 , x2 , , and xm input variables T 1 , T 2 , , and T m critical temperatures cy thermal yaw DX 2 and DX 4 translational errors in X-direction R distance between two sensors on the test rod DX 1 and DY 1 distances in X- and Y- directions from probe of temperature sensor S1 to test rod g Fitting ability/Predictive ability Tn temperature of moving nut cn and mn specific heat capacity and mass of nut Fa axial load A and P denote cross-section area and axial force, DL1 actual elongation DL thermal error with pre-stretch
five-axis CNC machines tools is reduced or avoided. So far, the thermally-induced error investigation becomes a hot topic of current research activities. The error avoidance method, the thermal characteristic optimization method, and the error compensation method are widely applied to eliminate the thermally-induced error. The error avoidance method attempts to eliminate or reduce possible error sources at the design and manufacturing stages [2]. The influence of error sources on the system is reduced with the increase in the design and manufacturing accuracy of machine tools’ components, and then the strict temperature control, the vibration isolation, the airflow disturbances, and the environmental state control are employed to eliminate or reduce the influence of external error sources of machine tools. Although it can reduce the original error, there still exist great limitations in the improvement of the manufacturing and installation accuracy of machine tools. It is shown that the cost of the error avoidance method increases exponentially when the machining accuracy requirements are higher than a certain level. Then the error avoidance method fails in improving the machining accuracy of machine tools costly. When the error avoidance method encounters difficulties, researchers switch their interest to the thermal and dynamic characteristic analysis and optimization of machine tools. The bearing thermal model [3–4], the thermal-structure coupling characteristics analysis [5–6], the thermal contact resistance (TCR) [7–9], the optimization of the fits for the bearing joints [10], the hydrodynamics analysis of machine tools [11], the contact stiffness of the tool/holder interface [12], and the precision loss modeling of machine tools [13] are involved in the thermal characteristic analysis of machine tools. Moreover, the thermal characteristic analysis is not only conducted for the whole machine tools, but also is conducted for critical functional components of machine tools [14–20]. The thermal characteristic analysis method is effective to analyze the temperature rise and the thermal deformation of machine tools. But accurate heat source intensity and thermal boundary conditions
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
5
are difficult to be obtained because of the nonlinear temperature rise and thermal displacements of machine tools, resulting in the modeling error of the thermal characteristic. Moreover, the results obtained by the thermal characteristic analysis are difficult to be used in the subsequent error compensation because the numerical modeling is time-consuming. The above factors limit the wide application of the thermal characteristic analysis and optimization methods. The thermally-induced error compensation is then becoming more and more popular in the past few decades. Moreover, it is revealed that the thermally-induced error is reduced significantly when the compensation method is implemented [21– 23]. The main process includes such steps as thermal error measurement, identification, modeling, and compensation. The thermal error measurement and identification are realized by the laser interferometer, etc. [24–25]. The thermal error model with excellent predictive performance and strong robustness is the core of the thermal error compensation according to previous studies [26–27]. Then various algorithms are used to establish the correlation between the thermal error and measured temperatures. The most widely-used algorithms in the thermal error modeling are the multivariate linear regression analysis (MLRA) [28], the time series (TS) [29], the genetic algorithm (GA) [30], and the artificial neural network (ANN) [31]. The error compensation of linear axes is also addressed to reduce the thermally-induced error [32–35]. However, the robustness of the above models needs further verification because of the dynamic change in the thermal behaviors, resulting from complex and volatile working conditions. The above models are the empirical error models with temperatures as the input and with the thermal error as the output, and the modeling method is less compatible with thermal error data without revealing the error generation mechanism. Moreover, for the five-axis machine tools, the error generation mechanism of has not been fully investigated to construct the thermal error model. Instead, the thermally-induced error is correlated with critical temperatures empirically, which splits the relationship between the error generation mechanism and the modeling method. Furthermore, the thermal error model without physical significance is difficult to be used in practical application, including the thermal errors based on GA and ANN. Besides, the combined effects of thermal errors of the linear and rotational axes are considered simultaneously, leading to the failure in the significant improvement of the machining accuracy of five-axis machine tools. Then the integrated thermal error modeling and compensation methods with strong robustness and excellent predictive performance are necessary to be investigated. To realize the integrated thermal error modeling and compensation, the homogeneous transformation (HT) method is used [36–41]. The thermal characteristic experiments are conducted offline, and then the error terms are identified. The error terms are substituted into the transformation matrix and these error terms generally are the static errors at thermal equilibrium state in the above models, failing in the dynamic modeling and compensation. Moreover, the mapping relation between the thermally-induced error and critical temperatures is constructed empirically because it is believed that the thermal error is the direct reflection of the temperature rise. But the thermal error generation mechanism of the linear and spindle axes is somewhat complex because of the interaction among such factors as the structure parameters, the working conditions, the ambient temperature, the material properties, the thermal load intensity and position, and the cooling effect, leading to the deterioration of the robustness. The relationship between the temperature rise, the thermal expansion, the thermally-induced error, and the resultant positioning error of linear axes has not been explored. So there still exists shortcomings in the HT method, and it is difficult to realize accurate and robust error compensation. The main goal of this study is to reduce the thermally-induced error of five-axis machine tools. The remarkable difference from previous studies is the integration of the thermal errors of linear and rotational axes with the aid of the homogeneous transformation, and the thermal errors of the linear axes are identified with the product of the polynomial function with temperature as its independent variable and the polynomial function with position as its independent variable from the error generation mechanism of linear axes. The advantage of the data-driven thermal error compensation method is its excellent effectiveness and strong robustness at different working conditions. The most significant contribution is the realization of the thermally-induced error compensation by the compensation model with excellent predictive performance and strong robustness. The main principle of the proposed method is expressing the compensation components as the matrix of error terms of the linear and spindle axes by analyzing the error transmission chain of five-axis machine tools. In Section 2, the thermal error model was established by calculating the homogenous transformation matrixes of the workpiece chain and tool chain of the five axis CNC machine tools, and then the compensation component was obtained for each axis. In Section 3, the thermal error terms of the linear and spindle axes were identified by the experiments. In Section 4, the thermal error compensation was implemented for a five-axis CNC machine tool to verify the validity of the proposed method. In Section 5, the effects of the heat source and pre-stretch of the lead screw on the thermal error are discussed. In Section 6, the conclusions obtained from previous analysis were summarized.
2. Thermally-induced error modeling method of five-axis machine tool 2.1. Thermally-induced error transfer model A standard five-axis CNC machine tool is made up of two rotational axes of B and C and three translational axes of X, Y, and Z, as depicted in Fig. 1. In the machining process, the cutter is installed on the spindle system in Z-direction, and the workpiece is installed on C-axis. The settings of the coordinate system are listed in Table 1 to assess the thermallyinduced error of the five-axis CNC machine tool. Three translational axes X, Y, and Z are perpendicular with each other mutually. The directions of the local coordinate systems are the same as that of the five-axis machine tool’s coordinate system. The
6
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 1. Configuration schematic view of five-axis CNC machine tools.
Table 1 Definition of coordinate systems. Coordinate system
Physical significance
x y z b
Local Local Local Local
c
Local coordinate system fixed to C-axis
s
Local coordinate system fixed to spindle system Coordinate system of tool nose Coordinate system of workpiece
t w
coordinate coordinate coordinate coordinate
system system system system
Directions fixed fixed fixed fixed
to to to to
X-axis Y-axis Z-axis B-axis
Remarks
Same as the direction of coordinate system of machine tools Coordinate origin O2 Coordinate origin O3 Coordinate origin O4
origins of the coordinate systems b, c, and s are O2 , O3 , and O4 , respectively. Moreover, the point O2 is on the B-axis, the point O3 is at the center point of the top surface of C-axis, and the point O4 is at the center point of the front end surface of the cutter. The reference coordinate system O1 is the intersection of B- and C- axes, and the intersection O1 is on the axis of S-axis. The transmission path of the thermally-induced error of five-axis machine tools is shown in Fig. 2. Moreover, the total transmission paths are divided into two transfer kinematic chains. One of the transmission paths of the kinematic chain is expressed successively as the machined workpiece!rotational axis C!rotational axis B!translational axis Y!translational axis X!machine tool bed. The other transmission path of the kinematic chains is expressed as the tool!the spindle system S!the translational axis Z!machine tool bed. The error will be transferred to each axis by the above two kinematic chains, resulting in the deterioration of the machining accuracy. The transformation matrixes between different coordinate systems are expressed as
7
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 2. Schematic view of kinematic chain.
2 b cT
cos a sin a 0 d2
6 sin a 6 ¼6 4 0
3
2
cos b
0
sin b
cos a
0
0
1
6 0 d1 7 1 0 7y 6 7 bT ¼ 6 4 sin b 0 cos b 0 5
0
0
1
0
0
0
d2
3
2
1 0
0
0
3
60 1 0 y7 07 7x 7 6 7 yT ¼ 6 7 40 0 1 05 05
0
0 0
1
0
ð1Þ
1
where a denotes the rotational angle of C-axis; b denotes the rotational angle of B-axis; y denotes the moving distance of Y-axis; d1 denotes the offset of O2 in the x-direction of the reference system O1 ; and d2 denotes the offset of O3 in the y-direction of the reference system O1 . Then the actual homogeneous transformation matrix is corrected with the geometric and thermal errors of the Y-axis.
2 x 0 yT
eYX
1
6e 6 YX ¼6 4 0
0
32
1 0
0
0
32
0 0
0 1
0
0
0
0 1
0 eyy dxy 1 0 0 6 0 1 0 Eðy; t; T Þ 7 exy dyy 7 7 76
1
exy
1
0
0
0
3
32
ezy
1
76 6 07 76 0 1 0 y 76 ezy 76 76 4 5 0 0 1 0 54 eyy 1 0
1
0
0
76 dzy 54 0 0
1
0 1
0
0
1
0
ð2Þ
7 5
where eYX denotes the squareness error between Y- and X-axes; Eðy; t; T Þ denotes the thermally-induced error of Y-axis; and dxy ; dyy ; dzy ; exy ; eyy and ezy denote six position-dependent geometric errors of the Y-axis. The actual homogeneous transformation matrix is corrected with the geometric and thermal errors of the X-axis.
2 O1 0 x T
1 0 0
x
32
1
6 0 1 0 0 76 e 76 zx 6 ¼6 76 4 0 0 1 0 54 eyx 0
0
0 1
eyx dxx 1 0 0 Eðx; t; T Þ 7 6 0 exx dyx 7 7 76 0 1 0
1
0
3
32
ezx
exx
1
76 dzx 54 0
0
0
1
0
0 1
0
0
1
0
ð3Þ
7 5
where x denotes the moving distance of the X-axis; dxx ; dyx ; dzx ; exx ; eyx and ezx denote 6 position-dependent geometric errors of the X-axis; and Eðx; t; T Þ denotes the thermal error of the X-axis. The actual homogeneous transformation matrix is corrected with the geometric and thermal errors of the Z-axis.
2
1 6 0 6 O1 0 z T ¼ 6 4 eZX 0
0 1
eZY 0
32
32
3
32
1 ezz eyz dxz 1 0 0 0 eZX 0 1 0 0 0 76 7 6 6 0 eZY 0 7 1 exx dyz 7 7 76 0 1 0 0 76 ezz 76 0 1 0 76 7 76 76 1 0 54 0 0 1 z 54 eyz exz 1 dzz 54 0 0 1 Eðz; t; T Þ 5 0
1
0
0
0 1
0
0
0
1
0
0
0
ð4Þ
1
where z denotes the moving distance of Z-axis; dxz ; dyz ; dzz ; exz ; eyz and ezz denote 6 position dependent geometric errors of the Z-axis; eZX and eZY denote the squareness error between Z- and X- axes and that between the Z- and Y- axes, respectively; and Eðz; t; T Þ denotes the thermal error of the Z-axis.
8
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
According to the transfer relationship of the thermal and geometric errors between two kinematic chains, the following relationship is obtained. O1 0 z T
z sT
s tT
¼ Ox 1 T 0
x 0 yT
y bT
b cT
ð5Þ
c wT
In the ideal case, the transformation matrix of coordinate systems w to c is expressed as c wT
¼ ðOz 1 T 0
where ts T , sz T , and
0 z O1 T
2 t sT
¼
O1 0 x T
s tT Þ
z sT
x 0 yT
y bT
b cT
¼
t sT
s zT
0 z O1 T
O1 0 x T
denote the inverse matrix of the st T , zs T , and
1
0
0
0
O1 0 z T ,
x 0 yT
y bT
respectively. Then st T , sz T , and
0
0 z O1 T
are expressed as
3
7 s 1 6 60 1 0 07 ¼6 7 tT 40 0 1 L5 0
ð6Þ
b cT
ð7Þ
0 1
where L denotes the length of the cutter.
3 1 0 0 0 60 1 0 0 7 7 6 ¼6 7 4 0 0 1 d3 5 2
s zT
¼
z 1 sT
0
0
0
ð8Þ
1
where d3 denote the offset of O4 in z-direction of the reference system O1 . 0 z O1 T
¼ ðOz 1 T 0 Þ
1
ð9Þ
The transformation matrix of coordinate systems t to w is expressed as w t T
O1 0 O1 0 z s w c b y 0x ¼E¼w O1 T t T ¼ c T b T y T x T O1 T z T s T t T
ð10Þ
where E denotes the identity matrix. 2.2. Thermally-Induced error modeling principles The thermal error terms of the spindle system include the angular errors of cx , cy , and cz and the translational errors of gx , gy , and gz , as shown in Fig. 3. Then the transformation matrixes of the spindle system around three transitional axes are obtained.
Fig. 3. Schematic graph of thermal error components of spindle system.
9
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
2
1
0
0
6 0 cosc 6 x Rð X; cx Þ ¼ 6 4 0 sincx 0
coscx 0
0 coscz
6 sinc 6 z RðZ; cz Þ ¼ 6 4 0
ð11Þ
1
0
sincy
1
0
0
3
0 coscy
07 7 7 05
0
1
0
sincz
0
0
ð12Þ
3
0
07 7 7 1 05
0
0 1
coscz
0
3
07 7 7 05
sincx
0
2 coscy 6 0 6 R Y; cy ¼ 6 4 sincy 2
0
0
ð13Þ
The spindle system rotates continuously around three coordinate axes of a fixed coordinate system for three times, and then
2
coscy coscz
6 sinc sinc cosc þ cosc sinc 6 x z z x z R ¼ Rð X; cx ÞR Y; cy RðZ; cz Þ ¼ 6 4 coscx sincy coscz þ sincx coscx 0
coscy sincz
sincy
sincx sincy sincz þ coscx coscz
sincx coscy
sincz sincz þ sincx coscz
coscx coscy
0
0
0
3
07 7 7 05 1 ð14Þ
Based on the assumption of small errors, then coscx ¼ coscy ¼ coscz ¼ 1, sincx ¼ cx , sincy ¼ cy and sincz ¼ cz . If the terms higher than the second-order are ignored, then Eq. (14) is simplified as
2
1
6 c 6 R ¼6 z 4 cy 0
cz 1
3
cy 0 cx 0 7 7
cx
1
7 05
0
0
1
ð15Þ
Combining Eq. (15) with the translational error of the spindle system, the transformation matrix movement of the tool nose coordinate systems t to w is expressed as
2 w t Te
1
6 c 6 ¼6 z 4 cy 0
cz 1
cx 0
twT e of the error
3
cy gx gx gy 7 7 7 1 gz 5 0
ð16Þ
1
Then the error compensation is achieved by the differential motion of each axis, and the error generated in the compensation movement is ignored. The compensation movement is the inverse movement of the error movement. The coincidence of the cutter tool nose coordinate system t and the workpiece coordinates system w is realized when the compensation movement is carried out. Namely, w w t Tct Te
¼E
ð17Þ
where twT c denotes the homogeneous transformation matrix of the compensation movement of the coordinate systems t to w. By using the differential transformation principle, the total differential relationship between coordinate systems t and w is obtained.
dw t T ¼
@w @wT @wT @wT @wT t T Dx þ t Dy þ t Dz þ t Da þ t Db @x @y @z @a @b
ð18Þ
where Dx, Dy, Dz, Da and Db denote the differential movement of the compensation of the translational axes of X, Y, and Z and the rotational axes of C and B. Moreover, the partial differential is expressed as y @w c b y 0 dx T O1 0 z s t T ¼w T T T c T bT y T x T @x dx z s t
ð19-aÞ
10
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538 y
dO1 T @w c b t T ¼w c T bT y T @y dx
0
0 O1 z s x O1 T z T s T t T
ð19-bÞ
@w dOz 1 T z s c b x y t T ¼w T T c T b T x T y T O1 T @z dz s t
ð19-cÞ
dcb T b y x O1 z s @w t T ¼w T T T T T T c T @a da y x O1 z s t
ð19-cÞ
dby T y x O z s @w c t T ¼w T T 1T T T c T bT @b db x O1 z s t
ð19-dÞ
The homogeneous transformation matrix of the compensation movement between the coordinate systems t and w is expressed as
2
1 6 k z 6 w t T c ¼6 4 ky 0
kz
ky
1 kx
kx 1
0
0
rx 3 ry 7 7 7 rz 5
ð20Þ
1
where kx , ky , and kz denote the angular errors in the compensation movement process, and rx , ry , and rz denote the transw lational errors in the compensation movement process. It is revealed that w t T c has an identical form with t T e . According to the error compensation principle, the following differential relationship is obtained. w t Tc
w ¼w t T þ dt T
ð21Þ
where dw t T denotes the change in the coordinate system caused by the differential transformation. Then the following relationship is obtained according to Eq. (21).
kx ¼ cx ; ky ¼ cy ; kz ¼ cz
ð22-aÞ
rx ¼ gx ; ry ¼ gy ; rz ¼ gz
ð22-bÞ
It is shown that the unique solution cannot be obtained by Eq. (22). The angular error cz of the tool cutter has no effect on machining accuracy of CNC machine tools because the rotational movement of the tool cutter relative to the spindle coordinate system z is symmetrical. Therefore, it is not necessary to compensate the thermal error in Z-direction of the workpiece coordinate system w and the cutter tool nose coordinate system t. Then the thermal error compensation components are obtained as
2 @ w T 31 t 3 2 3 2 3 gx rx Dx 7 6 @x @w t T 7 7 6 g 7 6 6 6r 7 7 6 @y 6 Dy 7 6 6 y7 y7 7 6 7 6 @w T 7 7 6 7 6 6 Dz 7 ¼ 6 gz 7 ¼ 6 6r 7 t 7 7 6 7 6 6 @z 7 6 z 7 6 7 6 7 6 6 7 w 7 4 Da 5 4 cx 5 6 6 @@t aT 7 4 kx 5 5 4 cy ky Db @w T 2
ð23Þ
t
@b
It is seen that the unique solution is obtained by Eq. (23). Then the relationship between the error compensation components Dx, Dy, Dz, Da, and Db are obtained, and then the thermal errors gx , gy , gz , cx , and cy of the spindle system are then obtained. Moreover, it is noticed that
@w T @w T @w T @w T t t t , @y , @z , @t a , @x
and
@w T t @beta
are closely related to geometric errors and thermal errors
of the linear axes. Then the following error terms are to be identified, including the geometric error terms dxx ; dyx ; dzx ; exx ; eyx , ezx and the thermal error term Eðx; t; T Þ of X-axis, the geometric error terms dxy ; dyy ; dzy ; exy ; eyy , ezy , and eYX , and the thermal error term Eðy; t; T Þ of Y-axis, the geometric error terms dxz ; dyz ; dzz ; exz ; eyz , ezz , eZX , eZY and the thermal error term Eðz; t; T Þ of Z-axis, and the thermal error terms cx , cy , cz , gx , gy and gz of the spindle system. The elimination of thermal error terms of the linear and spindle axes is the focus, including Eðx; t; T Þ, Eðy; t; T Þ, Eðz; t; T Þ, cx , cy , cz , gx , gy and gz . Then the identification of the geometric errors of exx ; eyx ; ezx , exy ; eyy ; ezy , exz ; eyz ; ezz , eYX , eZX , and eZY will be not studied anymore in the following section, and the identification of the thermal error terms of the linear and spindle axes, including Eðx; t; T Þ, Eðy; t; T Þ, Eðz; t; T Þ, cx , cy , cz , gx , gy , and gz will be investigated. As is shown, the compensation components in different directions are closely related to the position, time, temperature, the offsets d1 , d2 , and d3 , and the length L of the cutter. Moreover, the compensation component of each axis can be obtained according to Eq. (23) if the thermal error terms of the linear and spindle axes are identified.
11
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
2.3. Thermally-Induced error modeling method 2.3.1. Thermally-Induced error modeling of linear axis The thermal error modeling method of the linear axis will be studied from the view of the thermal error generation mechanism. The temperature field of the screw shaft is analyzed, and then the thermal error model of the linear axis will be proposed based on the thermal field modeling. (1) Heat Generation of Screw Nut Pair The heat generated by the motion of the screw nut pair is expressed as
Q2 ¼
2p n2 M2 60
ð24Þ
where M 2 denotes the friction torque between the nut and screw shaft; and n2 denotes the rotational speed of the screw shaft and can be expressed as
n2 ¼
V Ph
ð25Þ
where P h denotes the lead of the screw shaft; and V denotes the feed rate. The friction torque of the screw shaft pair is expressed as
M ¼ 2zðMe þ Mg Þcosb0
ð26Þ
0
where b denotes the helix angle of the raceway; and Me and Mg denote the friction torques caused by the external load and the geometric slip of balls, respectively. According to Ref. [20], Me and Mg are expressed as
sffiffiffiffiffiffiffiffiffiffiffiffiffi 4Q 4iq 3 P Me ¼ mb # q
ð27aÞ
vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi u 16Q 5iq 3 f m2a u t Mg ¼ 0:08 P R ð# qÞ2
ð27bÞ
where ma and mb denote the coefficients related to the eccentricity of the contact deformation ellipse, respectively; f denotes the coefficient of sliding friction; # denotes the parameter related to the material properties and is expressed as P q denotes the # ¼ 28 2 , wherein l and E denote the Poisson’s ratio and modulus of elasticity, respectively; and 3
1l 1l 1þ 2 E1 E2
sum principle curvatures of the contacting bodies at the contact point on the orthogonal plane and can be expressed as
X
q ¼ q11 þ q12 þ q21 þ q22
ð28Þ
a For the groove surface of the screw side, q11 ¼ q12 ¼ D2 , q21 ¼ 1R and q22 ¼ 2cosbcos , in which R denotes the raceway radius. dDcosb a For the groove surface of the nut side, q11 ¼ q12 ¼ D2 , q21 ¼ 1R and q22 ¼ 2cosbcos . dDcosb
The load imposed on the balls is shown in Fig. 4(a). As is shown, the ball is under the combined effect of the gyroscopic moment M gy , centrifugal force F cq , the shear loads F iq and F eq and the normal loads Q iq and Q eq . When the rotational speed of
Fig. 4. Quasi-static mechanics analysis of ball screw nut.
12
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
screw shaft is relatively high, the positions of the balls and the curvature centers of inner and outer groove are shown in Fig. 4(b). A1 and A2 denote the x and z coordinates of the inner and outer groove curvature centers and can be obtained according to Ref. [5].
(
A1 ¼ BDsina0 þ da þ hRi coswq
ð29Þ
A2 ¼ BDcosa0 þ dr coswq
where B denotes the inner and outer groove curvature constant; a0 denotes the initial contact angle; D denotes the diameter of the balls; h denotes the angular displacement of the inner ring; wq is the position angle of the qth ball; Ri denotes the radius of the groove curvature center; and da and dr denote contact deformations of the balls and grooves and are obtained according to the Pythagorean theorem and the geometric relationship shown in Fig. 4(b).
1=2 da ¼ V 21 þ V 22 ðf i 0:5ÞD h i1=2 ðf o 0:5ÞD dr ¼ ðA1 V 1 Þ2 þ ðA2 V 2 Þ2
ð30Þ
where f i and f o denote the curvature radius coefficients of inner and outer grooves, respectively; and V 1 and V 2 denote the distance variation between the curvature center of the inner ring and balls, respectively. The relationship between the center of the qth ball and the curvature centers of inner and outer grooves is expressed as
(
ðA1 X 1 Þ2 þ ðA2 X 2 Þ2 ¼ ½ðf i 0:5ÞD þ da X 21
þ
X 22
¼ ½ðf o 0:5ÞD þ dr
2
2
ð31Þ
where X 1 and X 2 denote the initial values of the x and z coordinates of balls, respectively. The qth ball is under the combined effect of the gyroscopic moment M gy , centrifugal force F cq , the shear loads F iq and F eq and the normal loads Q iq and Q eq . Then, the force balance equations of the qth ball is expressed as M gq cos eq d M gq sin eq iq d
Q oq sinaeq Q iq sinaiq Q oj cosaeq Q iq cosa
a cosaiq ¼ 0 a sinaiq þ F cq ¼ 0
ð32Þ
3=2 where Q iq ¼ K i d3=2 iq ; Q oq ¼ K o doq , wherein K i and K o denote the contact stiffness of the ball / nut raceway and ball / screw
raceway; aeq and aiq denote the contact angles of the balls / nut groove and balls / screw shaft groove, respectively. The principal curvature function is expressed as
F ðqÞ ¼
j2 þ 1 R 2C ðj2 þ 1ÞR
ð33Þ
P 0:636 P P P q q q where j ¼ 1:0339 P q1 ; C ¼ 1:5277 þ 0:6023ln P q1 ; and R ¼ 1:0003 þ 0:5968ln P q1 , wherein q1 ¼ q11 þ q22 2 2 2 P and q2 ¼ q12 þ q21 . The axial load F a , radial load F r , and moment load M g of the inner ring are expressed as
8 z P > > Fa Q iq sinaiq ¼ 0 > > > q¼1 > > > < z P Q iq cosaiq ¼ 0 Fr > q¼1 > > > > z P > > > Q iq Ri sinaiq ¼ 0 : Mg
ð34Þ
q¼1
where z denotes the number of balls. The model of the screw nut pair is W4009SA-6D-C5Z12 and is manufactured by the NSK Company. The basic parameters of the screw nut pair are as follows: the total length is 1269 mm, the diameter is 40 mm, the length of the nut is 225 mm, the maximum stroke is 700 mm and the lead is 12 mm. The rotational speeds of the screw shaft are 416.7 r/min, 666.7 r/min, 1000 r/min and 1250 r/min when the feed rates are 5 m/min, 8 m/min, 12 m/min and 15 m/min, respectively. The heat generation of the screw nut pair under different feed rates is obtained, seen in Fig. 5. As is shown, the heat generation of the screw nut pair increases with the feed rate. (2) Heat Generation of Bearings According to Ref. [6], the heat generated by the bearings is expressed as
(
Q 1 ¼ 1:047 104 n1 M1 M1 ¼ Ml þ Mv
ð35Þ
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
13
Fig. 5. Heat generation of screw nut pair.
where M 1 denotes the friction torque of the bearings; n1 denotes the rotational speed of the screw shaft; and M l and M v denote the friction torques caused by the external loads and the lubricant viscous friction, respectively. The friction torque M l caused by the external load is expressed as
M l ¼ f 1 F b dm
ð36-aÞ
F b ¼ maxð0:9F a =tana 0:1F r ; F r Þ
ð36-bÞ
where f 1 denotes the coefficient related to the bearing structure; F b denotes the equivalent load applied on bearings; dm denotes the pitch diameter of bearings; Fa and Fr denote bearing’s axial and radial loads, respectively; max denotes the greater value between 0:9F a =tana 0:1F r and F r ; and a denotes the contact angle of the bearing. The friction torque M v caused by lubricant viscous friction is expressed as
(
Mv ¼ 107 f o ðv o nÞ2=3 dm Mv ¼ 160 10
7
f o dm
3
3
v o n P 2000 v o n < 2000
ð37Þ
where v o denotes the lubricant kinematic viscosity at operating temperature; n denotes the rotational speed of the bearing; and f o denotes the constant related to the bearing model and lubrication type. f o is taken as 1 for single oil-mist or oil-gas lubrication. f o is taken as 2 for oil lubrication and is taken as 4 for grease lubrication of double column. The quasi-static analysis method, which is similar with the method shown in Ref. [3] to Ref. [6], is used to compute the heat generation and the friction torques of the bearings and the ball screw. The nonlinear balance equation set of bearings, including the balance equations of the geometrical, physical, rolling element, and bearing rings, are established. Then the Newton-Rapson algorithm is used to solve nonlinear balance equation set and the heat generation and the friction torques of the bearings are obtained. The parameters of the bearings are as follows. The inner diameter, the outer diameter, the thickness, and the contact angle, and the curvature radius coefficient of the grooves are 40 mm, 90 mm, 20 mm, 60°, and 0.55, respectively. Then the heat generation of bearings is obtained, seen in Fig. 6. It is depicted that the heat generation increases with the feed rate and the preload. (3) Convective coefficient According to Ref. [42], the convective coefficient is expressed as
14
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 6. Heat generation of bearings.
h¼
Nu k L
ð38Þ
where L denotes the feature size; k denotes thermal conductivity of fluid; and Nu denotes the Nusselt number. The heat transfer mode of moving components, including the screw shaft, belongs to forced convection. The criterion equation for the forced convection is expressed as
Nu ¼ 0:133Re2=3 Pr1=3
ð39Þ
where Re and Pr denote the Reynolds number and Prandtl number, respectively. The Reynolds number Re is expressed as
Re ¼ wd0 =v 2
ð40Þ
where w, v , and d0 denote the angular velocity, the viscosity of air, and the diameter of the screw shaft, respectively. The angular velocity is expressed as
w ¼ 2pn2
ð41Þ
where n1 denotes the rotational speed of the screw shaft. The forced convective coefficients of the screw shaft under different feed rates are obtained according to the above analysis, as shown in Table 2. It is seen that the forced convective coefficient increases with the fee rate. (4) Thermal Field Modeling The cross-section area and the perimeter of the shaft along the longitudinal direction do not change by assuming that the screw shaft is an extension rod with a circular cross-section. The temperature of the screw shaft at the axial position x and time t is T ðx; tÞ. Then the temperature rise of the screw at the different axial positions is expressed as T ðx; tÞ T 0 . According to Fig. 7, the heat transfer equation is obtained as
Table 2 Forced convective coefficients of screw shaft. Feed rate (m/min)
Rotational speed (r/min)
Forced convective coefficients (W= m2 C )
5 8 12 15
416.7 666.7 1000 1250
70.1 75.6 81.2 86.4
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
15
Fig. 7. Heat transfer of the shaft.
@ 2 T ðx; tÞ qc @T ðx; tÞ 4h þ dx ¼ ðT ðx; t Þ T 0 Þ @x2 k @t kd0
ð42Þ
where T 0 denotes the temperature of the ambient; h denotes the convective coefficient; and a0 denotes the thermal expansion coefficient; q and c denote the density and the specific heat capacity, respectively; d0 and T 0 denote the diameter of the screw shaft and the reference temperature, respectively; k denotes the thermal conductivity; x and t denote the position and the time, respectively. When the linear axis reaches the thermal equilibrium state, Eq. (42) is written as 2
d T ð xÞ 4h ¼ ð T ð xÞ T 0 Þ dx2 ad0 Assuming that
4h
ad0
ð43Þ
¼ m2 , then Eq. (43) is expressed as
2
d T ð xÞ ¼ m2 T ðxÞ dx2
ð44Þ
If the thermal expansion coefficient a is a constant, Eq. (42) is a quadratic linear differential equation with constant coefficients. Then the solution of Eq. (44) is expressed as
T ðx; t Þ ¼ C 1 emx þ C 2 emx
ð45Þ
Moreover, the following relationships are satisfied by Eq. (45). (1) T ¼ T ð0Þ when x ¼ 0, wherein T ð0Þ denotes the temperature of the heat source. ¼ aL T L , wherein L and T L denote the length and the temperature of the screw shaft, respectively. But T L is (2) k dT dx x¼L unknown, then it is necessary to add another condition to determine C 1 and C 2 . Namely, (3) T ¼ T L when x ¼ L. Then the following relationship is obtained.
T ð0Þ ¼ C 1 þ C 2
ð46Þ
C 1 memL C 2 memL ¼ aT L
Then C 1 and C 2 are obtained according to Eq. (46), and then the temperature distribution in the axial-direction is expressed as
T ðx; t Þ ¼
ðkm þ hÞemðLxÞ þ ðkm hÞemðLxÞ T ð0Þ ðkm þ hÞemL þ ðkm hÞemL
ð47Þ
The temperature rise stage is the none-equilibrium state. The efficiency of the heat conduction is far higher than that of the convection. Then the heat transfer equation is simplified as
@ 2 T ðx; tÞ qc @T ðx; tÞ dx ¼ @x2 a @t
ð48Þ
By solving Eq. (48), the temperature response is then obtained.
"
2 T ðx; t Þ ¼ T ð0Þ 1 pffiffiffiffi
p
Z px ffiffiffiffi a 2
qct
#
e
u2
du
ð49Þ
0
Then the temperature model of the screw shaft under the effect of a single heat source is obtained according to Eq. (45) and Eq. (49).
16
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
T ðx; tÞ ¼ T ð0Þ
ðkm þ hÞemðLxÞ þ ðkm hÞemðLxÞ 2 1 pffiffiffiffi ðkm þ hÞemL þ ðkm hÞemL p
Z px ffiffiffiffi a 2
qct
! 2
eu du
þ T 0 ¼ T ð0Þ T ðxÞ T ðt Þ þ T 0
ð50Þ
0
where T ðxÞ and T ðt Þ denote the position-dependent and time-dependent components of the temperature, respectively. It is necessary to know the parameters in Eq. (50) to get the thermal field. Then the parameters in Table 3 are used to compute the thermal field at the thermal equilibrium state, seen in Fig. 8. It is shown that the temperature increases with time and that the temperature does not change significantly with position at the temperature equilibrium state. However, the screw shaft is under the combined effect of three heat sources, including the front bearing, the rear bearing, and the nut. Then Eq. (50) is not the exact solution of the temperature response of the screw shaft. Instead, the exact solution is the superstition of the temperature responses of three heat sources. Namely,
T ðx; tÞ ¼ T f ðx; tÞ þ T r ðx; tÞ þ T s ðx; t Þ þ T 0
ð51Þ
(5) Thermal Error Modeling and Identification Methods For the linear axis with the fixed-support or fixed-free installation, the transient thermal elongation, caused by the temperature rise, is then obtained as
RL ða0 DT ðx; t ÞÞdx ¼ 0 ða0 ½T ðx; t Þ T 0 Þdx 0 RL 0 ¼ 0 a T f ð0Þ T f ðxÞ T f ðt Þ þ T r ð0Þ T r ðxÞ T r ðt Þ þ T s ð0Þ T s ðxÞ T s ðtÞ dx
DLðt; x; T Þ ¼
RL
ð52Þ
where DT ðx; t Þ denotes the temperature rise; a0 denotes the thermal expansion coefficient; T f ð0Þ, T f ðxÞ, and T f ðt Þ denote the initial temperature, the position component, and the time component of front bearing’s temperature; T r ð0Þ, T r ðxÞ and T r ðtÞ denote the initial temperature, the position-dependent component, and the time-dependent component of rear bearing’s temperature; and T s ð0Þ, T s ðxÞ and T s ðtÞ denote the initial temperature, the position-dependent component, and the timedependent component of nut’s temperature. The Taylor expansion is conducted for both the temperature components T ðxÞ and T ðtÞ, and then they are obtained as
T ðxÞ ¼
1 X T ðkÞ ð0Þ k x k! k¼0
ð53-aÞ
T ðt Þ ¼
1 X T ðkÞ ð0Þ k t k! k¼0
ð53-aÞ
By substituting Eq. (53) into Eq. (52), it is revealed that the thermal error is approximated by the polynomials with time and position as their independent variables. Then the thermal error of the linear axis is approximated by these polynomials. It is noticed that the thermal expansion at different positions is a function with time t as its independent variable. Then it is difficult to calculate the thermal expansion of the screw shaft directly by Eq. (52) and Eq. (53). But the above analysis provides the basic idea for the thermal error fitting of the linear axis. If the high-order terms are ignored, then the thermal error of the linear axis is expressed as
Eðx; t; T Þ ¼ DLðt; x; T Þ ¼ At t 3 þ Bt t þ C t ax x3 þ bx x þ cx
ð54Þ
where At , Bt , and C t , which are related to the time, are the coefficients to be determined; and ax , bx , and cx , which are related to the position, are the coefficients to be determined. The temperature of the screw shaft is regarded as the intermediate variables to simplify the decomposition of the thermally-induced error from the measured positioning error of linear axes. The measured positioning errors of linear axes are the sum of the thermal error and the geometric error, namely,
DEðx; t; T Þ ¼ Eðx; t; T Þ þ oðxÞ
ð55Þ
where oðxÞ denotes the geometric error. The thermal error Eðx; t; T Þ is a function with time and position as its independent variables, while the geometric error oðxÞ is a function with position as its independent variable. The above operation greatly simplifies the thermal error modeling of the linear axis, and then the decomposition of the thermal error from the measured positioning error is realized. The basic idea for the modeling method of thermal error of the linear axis is then obtained.
Table 3 Parameters used in Eq. (50). Parameters
a (m/°C)
h (W/m2/ °C)
k (W/(m °C))
T 1 (°C)
T ð0Þ (°C)
L(mm)
Values
1.0e-5
70.1
60
21
60
400
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
17
Fig. 8. Temperature distribution of the linear axis.
2.3.2. Thermally-induced error modeling of spindle system The thermal error of the spindle system is closely related to temperatures of critical points because the thermal error is the manifestation of the temperature rise. The MLRA model is widely adopted because of its simple structure and highefficiency calculation. Moreover, the MLRA-based thermally-induced error model has a clear physical meaning, and then it is extremely convenient for engineering applications. The thermal error modeling method of the spindle system will be studied based on the MLRA. (1) Thermally-induced Error Modeling Based on MLRA Method The MLRA is effective to establish the correlation between the thermal error and the critical temperatures. Then the relationship between the dependent variable and the independent variable is
y ¼ b 0 þ b 1 x1 þ b 2 x2 þ þ b m xm þ e
ð56-aÞ
where b0 , b1 , b2 , bm denote the coefficients to be determined; e denotes the residue; x1 , x2 , , and xm denote the input variables; and y denotes the output variable. The least-square method is used to calculate the coefficients. The thermal error models of the spindle unit, including the thermal elongation Ez , thermal yaw cx , thermal pitch cy , the translational error gx and gy , are expressed as
0
Ez
1
0
b0E
Bc C Bb B x C B 0cx B C B B cy C ¼ B b0cy B C B Bg C B @ x A @ b0gx gy b0gy 0
b1E
b2E
b1cx
b2cx
b1cy
b2cy
b1gx
b2gx
b1gy
b1gy
10 1 1 0 1 eEz C B T1 C B bmcx C ec C CB C C B CB B xC C B T C B bmcy CB 2 C þ B ecy C C CB .. C C C B bmgx AB @ . A @ egx A egy b1gy Tm bmE
ð56-bÞ
1 0 1 b0E b1E b2E bmE eEz B b0cx b1cx b2cx bmcx C B ecx C B C B C C B C where B B b0cy b1cy b2cy bmcy C denotes the coefficient matrix; B ecy C denotes the residual matrix; and T 1 , T 2 , , and @ b0g b1g b2g bmg A @ egx A x x x x egy b0gy b1gy b1gy b1gy T m denote the critical temperatures. (2) Angular errors According to Fig. 9, the translational error of the spindle system in the Z-direction is measured by the displacement sensor S5. Moreover, the thermal yaw cy is measured by displacement sensors S1 and S3. Namely,
18
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 9. The spindle thermal inclination sketch.
tancy ¼
DY 1 DY 3 DH
ð57-aÞ
where DY 1 and DY 3 denote the translational errors of the spindle system in Y-direction; and DH ¼ 100mm. Then the thermal yaw cy is expressed as
DY 1 DY 3 DH
cx ¼ arctan
ð57-bÞ
Correspondingly, the thermal pitch cx is expressed as
tancx ¼
DX 2 DX 4 DH
ð58-aÞ
where DX 2 andDX 4 denote the translational errors of the spindle system in X-direction. Then, the thermal pitch cx is expressed as
DX 2 DX 4 DH
cx ¼ arctan
ð58-bÞ
(3) Translational errors The translational error of the spindle system in the Y-direction is mainly attributed to the angular error around Zdirection and that in Y-direction. So the translational error in Y-direction is expressed as
DX 1 ¼ gx þ R cos/ ð1 coscy Þ
ð59-aÞ
where / and R denote the rotational angle of the spindle system and the distance between two sensors on the test rod, respectively. Similarly,
DY 1 ¼ gy þ R cos/ ð1 coscx Þ
ð59-bÞ
where gy denotes the translational error in the Y-direction. So the translational errors gx and gy of the spindle system in the X- and Y-directions are then identified according to Eq. (59).
(
gx ¼ DX 1 R cos/ ð1 coscy Þ gy ¼ DY 1 R cos/ ð1 coscx Þ
ð60Þ
where DX 1 andDY 1 denote the distances in X- and Y- directions from the probe of temperature sensor S1 to the test rod, respectively.
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
19
3. Thermally-Induced error identification of linear and spindle axes 3.1. Experimental setup and measuring principles The thermal error terms, including the Eðx; t; T Þ, Eðy; t; T Þ, Eðz; t; T Þ, gx , gy , gz , cx , cy , and cz will be accurately characterized. The angular errors and the translational errors of the spindle system are measured. The angular errors of the spindle system were measured by the five-point method [43]. The displacement senor S5 was used to collect the axial thermal elongation. The displacement sensors (S1 , S3 , S2 and S4 ) were used to measure the thermal yaw cx and thermal pitch cy , as shown in Fig. 10. Then the National Instrument cRIO-9045 was used to construct the acquisition system, in which high-precision C8-2.0 displacement and PT100 sensors were used to measure the thermal errors and the temperature rise, respectively. The performance parameters of the displacement sensors are listed in Table 4. The temperature sensors are used to capture the temperatures of the main heat sources of the spindle system. The temperature sensors, which are the OMEGA PR-11-2-100-M15-(*)-E, are the PT100 sensors. The temperature sensor is installed on the principle of capturing the temperatures of the main heat sources, and then the installation positions of these temperature sensors are listed in Table 5. For the spindle system, the thermal error is closely related to the critical temperatures T 1 ~T 5 . The laser interferometer Renishaw XL80 was used to collect the positioning errors of the X-, Y-, and Z- axes. Then the thermal error will be decomposed from the measured positioning error. The positioning error of the linear axis with fixed-support installation was investigated, depicted in Fig. 11. To simulate the thermal equilibrium process, the nut with the worktable travels back and forth in the full stroke of the linear axis. Then the Renishaw laser interferometer XL80 was applied to collect the positioning error at each measurement point. The measuring principle is shown in Fig. 12. The whole stroke of the linear axis is from 0 ~ 700 mm. Then the measuring range of the Renishaw laser interferometer XL80 is set as 0~650 mm because of the installation of the laser interferometer. Here the reference point is set as X = 0 and the measuring spacing is set as 50 mm. So the total number of the measuring points is 11. For each measurement, the Renishaw laser interferometer XL80 takes about 2 s to accomplish the measurement and collection. Then the worktable with the retroreflector stagnates for 4 s to guarantee the complete the above process. The reverse overrun was set as 2 mm to avoid the measurement error caused by the backlash. The linear axis continuous operating for 9 h to reach the thermal equilibrium state and then the worktable moves back and forth for several times in the full stroke of the linear axis. The measured positioning errors are obtained as the average value of the collected data. The air pressure and humidity sensors were used to compensate for the fluctuation effect of the air pressure and humidity. The feed rate, the rotational speed, and the ambient temperature are the main factors, which have a significant effect on the thermal error of the five-axis machine tool. By changing the above three factors, then the working conditions of the fiveaxis machine tool are set according to Table 6 to simulate the actual machining process. The working conditions #1, #2, and #3 are used for the modeling, the prediction and the verification, respectively.
Y-axis
X-axis Test bar S 3 S1 S5 S4
Spindle
S2
Worktable
Z-axis Fig. 10. Diagram of displacement sensor installation.
20
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Table 4 Performance parameters of displacement sensor. Type
Manufacturer
Model
Measuring range
Bandwidth
Linearity
Peak resolution
Capacitive
Lion
C8-2.0
250 mm
10KHz
0:15% FSO
6.32 nm
Table 5 Installation of temperature sensors. Temperature sensors
Installation positions
T1 T2 T3 T4 T5
Motor stator Front bearing housing Rear bearing housing Coolant outlet Spindle Housing
5-axis machine tool
Renishaw laser interferometer XL80 PC Fig. 11. Experimental setup.
Fig. 12. Measuring principle of positioning error of ball screw feed drive system.
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
21
Fig. 13. Rotational speed distribution.
3.2. Thermally-Induced error terms identification 3.2.1. Identification of Thermally-Induced error terms of linear axis It takes 9 h to conduct the thermal characteristic experiments of the linear axis, and the measured positioning errors of X-, Y-, and Z- linear axes are then obtained, seen in Fig. 14. It is depicted that the positioning error increases with the position in the negative direction and that the increase rate declines with the position. Moreover, the positioning errors at a different time have an identical shape, which means the measurement of the positioning error is reliable. Besides, the positioning error of the linear axis increases with time because the heat conduction efficiency is higher than that of the heat convection. It is revealed that the positioning error of the linear axis shows strong nonlinearity respective to the position. The empirical model is used in previous studies to establish the correlation between the thermal error and the temperature rise, which ignores the thermal error generation mechanism of the linear axis. Then the necessity to construct the thermally-induced error correlation of linear axes from the view of the thermal error generation mechanism becomes urgent according to Section 2.3.1. For the linear axis with fixed-support installation, the thermal error is formulated as the product of the polynomial function with time as its independent variable and the polynomial function with position as its independent variable. Then the thermally-induced error correlation is used to fit the measured data of X-, Y-, and Z-axes according to Eq. (54). Here the polynomial model with degree 3 is used, and the maximum fitting goodness of R2 is 0.95, which verifies the fitting performance of the proposed model in this paper. Then the fitted positioning errors are compared with the measured data, seen in Fig. 15. It is shown that the fitted positioning errors match well with the measured data in the whole moving range and the whole thermal equilibrium process, and then the validity of the proposed method is verified. The thermal error is the main component of the positioning error of the linear axis, and then the decomposition of the thermal error is feasible according to Section 2.3.1. The geometric error oðxÞ is the difference between the measured positioning error and the fitted thermal error, and the geometric error oðxÞ at t ¼ 0 is greater than that of at other time because the thermal effect has not come into effect at t ¼ 0. To exclude the errors in the measurement and modeling, the thermal error model is modified as
22
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Table 6 Working conditions of linear axis. Working conditions
#1 #2 #3
Parameters
Functions
Feed rate of linear axis
Rotation speed of spindle
Ambient temperature
24 m/min 30 m/min 18 m/min
Fig. 13(a) Fig. 13(b) Fig. 13(c)
22:4 0:5o C 19:2 0:3o C 19:2 0:5o C
Modeling Prediction Compensation
Fig. 14. Positioning errors linear axis.
ax x3 þ bx x þ cx Eðx; t; T Þ ¼ At t 3 þ Bt t þ C t a 0 x3 þ b 0 x þ c 0
ð61Þ
where a0 , b0 , and c0 , which are related to the thermal error at x ¼ 0, denote the coefficients at x ¼ 0. Finally, the correlation between the thermally-induced error and position is established, and the thermally-induced errors of X-, Y-, and Z-axes are identified according to the above analysis. Then the identified thermal errors are then used for subsequent compensation. Then the fitting goodness of the RMSE, the fitting ability g, the absolute mean, the absolute maximum, and the absolute minimum are used to evaluate the fitting performance of these models, seen in Tables 7–9. The fitting ability g for the X-axis at 0 h, 2 h, 5 h, 7 h, and 9 h are 94.99%, 93.27%, 94.36%, 93.71%, and 93.52%, respectively. The fitting ability g for the Y-axis at 0 h, 2 h, 5 h, 7 h, and 9 h are 96.13%, 96.11%, 98.00%, 98.67%, and 98.72%, respectively. The fitting ability g for the Z-axis at 0 h, 2 h, 5 h, 7 h, and 9 h are 97.87%, 97.19%, 98.20%, 98.64%, and 98.08%, respectively. The fitting performance of the thermal error models of the feed drive systems is excellent.
23
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 15. Comparison between fitted thermal error and measured positioning errors of linear axis.
Table 7 Model evaluation of fitting performance for X-axis. Indicators
RMSE Fitting ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time(h) 0
2
5
7
9
0.0689 94.99 0.2339 0.3686 0.0129
0.1201 93.27 0.3086 0.4843 0.0150
0.2462 94.36 0.4416 0.6883 0.0185
0.3710 93.71 0.5423 0.8496 0.0259
0.5366 93.52 0.6523 1.0244 0.0328
0
2
5
7
9
0.0200 96.13 0.1273 0.2000 0.0000
0.0827 96.11 0.2636 0.4000 0.0000
0.0618 98.00 0.2364 0.3999 0.1000
0.0418 98.67 0.2000 0.3000 0.1000
0.0491 98.72 0.2182 0.3000 0.2000
Table 8 Model evaluation of fitting performance for Y-axis. Indicators
RMSE Fitting ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time(h)
24
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Table 9 Model evaluation of fitting performance for Z-axis. Indicators
RMSE Fitting ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time (h) 0
2
5
7
9
0.0282 97.87 0.1555 0.2000 0.1000
0.0964 97.19 0.2909 0.5000 0.1000
0.0618 98.20 0.2405 0.3000 0.1000
0.0564 98.64 0.2182 0.4000 0.2000
0.1382 98.08 0.3636 0.5000 0.0282
Fig. 16. Temperature rise at working condition #1.
3.2.2. Identification of Thermally-Induced error terms of spindle system The critical temperatures of the spindle system at working condition #1 are obtained, as shown in Fig. 16. The trends of temperature rise and drop are consistent with the trend of the rotational speed distribution shown in Fig. 13 (a), which demonstrates the dependence of the temperature on the rotational speed. The heat generation of bearings increases with the rotational speed according to Ref. [6]. Then the temperature is the direct reflection of the heat generation of bearings. Finally, it is concluded that the temperature rise is closely related to the rotational speed distribution of working condition #1. Moreover, the fluctuation of the temperatures means that the change in speeds causes the difficulty to achieve a heat balance. Then the thermal deformations of the high-speed spindle at working condition #1 are obtained, as shown in Fig. 17(a). The thermal deformations increase with the rotational speed of the spindle system, and there exists a strong lag effect between the thermal deformation and the rotational speed. Then the thermal yaw cx and thermal pitch cy are obtained, as shown in Fig. 17(b). Finally, the translational errors gx and gy are identified according to Eq. (60), seen in Fig. 17(c). The critical temperatures shown in Fig. 17(a) are regarded as the inputs of the MLRA model, and then the thermal errors are predicted by the critical temperatures with the aid of the MLRA method, seen in Eq. (62). The residues between the measured data and fitting errors are in the range of [-0.49 mm, 0.5 mm]. For transitional errors, the fitting goodness of R2 is greater than 0.94, which verifies the effectiveness of the modeling method.
Ez ¼ 89:6619 þ 0:5752T 1 83:8363T 2 þ 10:6994T 3 þ 14:1400T 4 12:6225T 5
cx ¼ 31:1984 0:5500T 1 12:6472T 2 þ 7:8383T 3 þ 9:6489T 4 3:0353T 5 cy ¼ 78:1968 þ 2:6834T 1 þ 22:0344T 2 21:7061T 3 18:7970T 4 þ 12:4126T 5 gx ¼ 177:2860 0:5828T 1 11:2493T 2 þ 11:5722T 3 þ 4:3105T 4 þ 3:8117T 5 gy ¼ 72:7147 0:2376T 1 4:5866T 2 þ 4:7182T 3 þ 1:7575T 4 þ 1:5541T 5
ð62Þ
Then Table 10 lists the fitting goodness for the evaluation of the predictive performance. The fitting ability g of the models for Ez , cx , cy , gx , and gy are 96.28%, 96.10%, 96.44%, 94.47%, and 89.10%, respectively. The fitting performance of the thermal error models of the spindle system is excellent. 4. Modal validation 4.1. Thermally-Induced error prediction 4.1.1. Prediction of Thermally-Induced error terms of linear axis To verify the validity of the thermal error modeling method, the static thermal error models established in Section 3.2.1 are then used to predict the thermal errors of the linear axis at working condition #2, and the predicted thermal errors are
25
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 17. Thermal error of the spindle system at working condition #1.
Table 10 Model evaluation of fitting performance. Indicators
E (lm)
cx (00 )
cy (00 )
gx (lm)
gy (lm)
RMSE Fitting ability g (%) Mean of absolute Absolute maximum Absolute minimum
0.5815 96.28 0.7148 1.0417 0
0.5783 96.10 0.7114 1.0523 0
0.5802 96.44 0.7137 1.0382 0
2.9336 94.47 1.3662 4.5612 0.0204
4.1275 89.10 1.6662 4.3602 0.0274
compared with the measured data, as shown in Fig. 18. The residues between the predicted errors and the measured data are in the range of [-0.6 mm, 0.6 mm]. The predicted errors match well with the measured data in the whole measuring range, which validates the effectiveness and robustness of the thermal error modeling method. Then the prediction goodness of the RMSE, the fitting ability g, the absolute mean, the absolute maximum, and the absolute minimum are used to evaluate the fitting performance of these models, seen in Tables 11–13. The predictive ability g of these models has no significant decrease, which validates that the models are robust and generalizable. The excellent predictive performance and strong robustness come from the mechanism error modeling of the feed drive system.
4.1.2. Prediction of thermally-induced error terms of spindle system The thermal characteristic of the spindle system at working condition #2, including the temperature rise and the thermal deformations, was measured, seen in Fig. 19. It is depicted that both the temperature and thermal deformations show a strong nonlinear relationship with the running time. Moreover, the dependence of the temperature and thermal deformations on the rotational speed distribution is still significant. The maximum axial thermal elongation reaches about 50 mm, and the maximum radial deformation is about 40 mm.
26
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 18. Comparison between fitted thermal error and measured positioning errors of linear axis at working condition #2.
Table 11 Model evaluation of prediction performance for X-axis. Indicators
RMSE Predictive ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time(h) 0
2
5
7
9
0.4144 88.57 0.5649 1.0825 0.1760
0.5703 90.24 0.6287 1.2870 0.1014
0.5502 93.41 0.6250 1.2684 0.1000
0.3571 95.52 0.5370 1.1765 0.1512
0.7028 94.83 0.6592 2.0923 0.0127
The thermal errors of the spindle system are then predicted with the critical temperatures in Fig. 19(a) as the inputs. The predicted thermal errors and the measured data are then compared, seen in Fig. 20. It is revealed that the thermal error models established in Section 3.2.2 can predict the thermal errors of the spindle system accurately, which validates the effectiveness and robustness of the proposed method. Then Table 14 lists the predictive goodness for the evaluation of the predictive performance. The predictive ability g of the models for E, cx , cy , gx , and gy are 98.42%, 92.10%, 95.97%, 95.90%, and 96.27% respectively. The predictive performance of the thermal error models is great. Comparing with the fitting ability g, then the robustness and generalizable ability are validated.
27
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538 Table 12 Model evaluation of prediction performance for Y-axis. Indicators
RMSE Predictive ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time(h) 0
2
5
7
9
0.0355 94.41 0.1495 0.3366 0.0074
0.0269 96.74 0.1230 0.3837 0.0162
0.0236 97.86 0.1198 0.3215 0.0147
0.0367 98.28 0.1437 0.3432 0.0156
0.0291 98.83 0.1438 0.2722 0.0077
0
2
5
7
9
0.1053 92.01 0.2666 0.5190 0.0136
0.0509 96.22 0.1780 0.4500 0.0247
0.5094 91.02 0.6084 1.1661 0.1310
0.0363 98.66 0.1349 0.3473 0.0033
0.0634 98.49 0.2066 0.4478 0.0231
Table 13 Model evaluation of prediction performance for Z-axis. Indicators
RMSE Predictive ability g (%) Mean of absolute (lm) Absolute maximum (lm) Absolute minimum (lm)
Time(h)
Fig. 19. Thermal characteristics of spindle system at working condition #2.
4.2. Thermally-induced error compensation To improve the machining accuracy of five-axis machine tools, the thermal error compensation experiments were conducted on the CNC machine tools. The thermal error compensation principle of the spindle system is shown in Fig. 21. The thermal error compensation system includes the thermal characteristic acquisition and error compensation modules. The PT100 temperature sensors are used to collect the temperatures. The running time and actual position are collected from the CNC system. Then these signals are filtered, amplified, and A/D converted to get the actual physic information of machine tools. Then the thermally-induced error terms are identified. For the spindle system, the thermally-induced error is identified as a function with critical temperatures as its independent variables. For the linear axes, the thermal errors are identified as the product of the polynomial function with time as its independent variable and the polynomial function with position as its independent variable. Finally, the thermal error compensation components are calculated in each direction according to Eq. (23). The compensation components of the axis are sent to the CNC controller by the compensation module. The reverse of the error is superimposed to the output commands of the servo controller to offset the position deviation by comparing the actual and command positions of the worktable. Finally, the real-time thermal error compensation is realized. The five-axis machine tool is warmed up at working condition #3. Then thermal error compensation is carried out for linear axes. Moreover, the positioning errors of the linear axes are measured separately. The compensation results are then compared for X-, Y-, and Z- axes, as shown in Fig. 22. It is shown that the present error compensation can effectively reduce
28
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 20. Thermal error of spindle system at working condition #2.
Table 14 Model evaluation of prediction performance. Indicators
E (lm)
cx (00 )
cy (00 )
gx (lm)
gy (lm)
RMSE Predictive ability g (%) Mean of absolute Mean of absolute maximum Mean of absolute minimum
0.3234 98.42 0.4836 0.9937 0.0289
1.0561 92.10 0.9412 1.9351 0.1216
0.3205 95.97 0.4966 0.9441 0.0046
1.2083 95.90 0.9403 1.9589 0.0086
0.5946 96.27 0.6602 1.4273 0.0000
Fig. 21. Schematic diagram of compensation for thermal error of spindle system.
the positioning error, and the traditional pitch error compensation can also reduce the positioning error to some extent, which verifies the validity of the proposed method. The positioning error with the present error compensation is much smaller than that with the traditional pitch error compensation, which indicates that the necessity for the consideration of dynamic changes in the thermal error. For the present error compensation, the positioning errors of X-, Y-, and Z-axes are in the range of [4.8 lm, 3.5 lm], [5.5 lm, 3.0 lm], and [4.8 lm, 1.9 lm], respectively. For the traditional pitch error
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
29
Fig. 22. Comparison of error compensation at working condition #3.
compensation, the positioning errors of X-, Y-, and Z-axes are in the range of [9.8 lm, 9.7 lm], [10.9 lm, 8.1 lm], and [11.8 lm, 6.2 lm], respectively. The positioning errors of X-, Y-, and Z-axes are in the range of [4.8 lm, 3.5 lm], [5.5 lm, 3.0 lm], and [4.8 lm, 1.9 lm], respectively. For the linear axes with thermal error compensation, the positioning error is reduced significantly, which verifies the effectiveness of the thermal error compensation method and validates the statement that the thermal error is the main component of the positioning error further. It is concluded that the positioning errors of the linear axes are reduced obviously and that the thermally-induced error compensation method is effective and robust. Then the thermally-induced error compensation is implemented for the spindle system at working condition #3. The compensation results are then obtained, as shown in Fig. 23. It is shown that the thermally-induced errors of the spindle system are decreased significantly. The thermal elongation, the thermal yaw cx , and thermal pitch cy are in the range of [3.4 lm, 5.0 lm], [3.000 , 2.500 ], and [3.000 , 3.000 ], respectively. It is noticed that the thermal errors start to decreases from about 275 min, which is slightly lagging behind the decrease in the rotational speed shown in Fig. 13(c). 4.3. Actual Machining The validity of the data-driven thermal error compensation method is fully validated when actual machining is performed. So a large number of test samples were machined by the five-axis machine tool with and without thermallyinduced error compensation at working condition #3. Before the actual machining of test specimens, the motion simulation is conducted for the five-axis machine tool to simulate the actual machining process. The SS304 is chosen as the test material because of its wide use, as shown in Fig. 24. The workpiece is fixed on the B-axis and then the linkage of the X-, Y-, Z-, C-, B-, and S- axes realize the precision machining of workpiece. The feed rates of X-, Y-, and Z- axes are 1.2 m/min. The rotational speeds of the C-axis and B-axis are 0. The rotational speed of the spindle system is 5000 r/min. The cutting depth is 50 mm. It takes about one day to machine 32 test specimens. When the machining of all test specimens is finished, these specimens are
30
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 23. Comparison of the compensation results at working condition #3.
D2 D1
Fig. 24. Test samples.
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
31
Fig. 25. Machining error of the test samples.
placed in the laboratory with a constant temperature of 20 °C for 8 h, and then the roundness error of the test specimens is measured by the coordinate measuring machine, as shown in Fig. 25. The sizes of D1 and D2 are the most important parameters for the test specimens because of the assembly requirement. The basic dimensions of D1 and D2 are 70 mm and 95 mm, respectively. Then the roundness error of the workpiece is measured with the coordinate measuring machine. It is seen that the roundness error distribution is much closer to the zero line and that the roundness error is reduced by more than 85% and 37% with the present error compensation compared with that without compensation and that without traditional error compensation, respectively. Then the effectiveness of the thermally-induced compensation method is totally verified. Moreover, the trend of the roundness error without compensation is similar to the preload relaxation process caused by the temperature rise of the linear axis. 5. Results and discussion 5.1. Effect of heat source on thermal error Then the temperature of the moving nut and the bearings are obtained as
T n ¼ cnQm2 n þ T 0
ð63Þ
T f 1 ¼ T f 2 ¼ c Qm1 þ T 0 b
b
where T f 1 and T f 2 denote the temperatures of the bearing #1 and bearing #2, respectively; T n denotes the temperature of the moving nut; cb and mb denote the specific heat capacity and mass of the bearing, respectively; and cn and mn denote the specific heat capacity and mass of the nut, respectively. According to Eq. (24), Eq. (35), and Eq. (63), the temperatures of the moving nut and the bearings are obtained. Two ball screws, whose models are W4009SA-6D-C5Z12, are used in the dual-drive servo feed drive system. The material of the screw shaft is ASTM E52100. The density q, the specific heat capacity c, and the thermal expansion coefficient k are obtained. The diameter d0 of the screw shaft is 85 mm. Then the forced convective coefficient h is calculated according to Eq. (38), as listed in Table 2. Then the temperature of the screw shaft is calculated, seen in Fig. 26. The measured data are
Fig. 26. Temperature response of screw shaft.
32
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
Fig. 27. Positioning error of feed drive system.
compared with the fitting results. It is shown that the calculated temperature agree well with the measured data and that the temperature response of the screw shaft increases with the running time t, which is consistent with the results obtained by Ref. [44]. The temperature rise T ðx; tÞ of the screw shaft is obtained. Then the temperature rise T ðx; tÞ is substituted into Eq. (52) to calculate the thermal error of the feed drive system, seen in Fig. 27. It is seen that the calculated thermal error agrees well with the scatter points, which validates the effectiveness of the thermal error model. It is concluded that the effect of the heat source on the thermal error is realized by changing the temperatures of the heat sources. Then the heat entering linear axes is changed. By substituting the change in the heat generation into the thermal field and thermal error models of linear axes, the temperature and thermal error are predicted accurately. 5.2. Effect of pre-stretch of lead screw on thermal error When the ball screw feed system is working, a large amount of heat is generated, and the heat causes the thermal expansion of the screw shaft, and hence the lead is increased, leading to the decrease in the positioning accuracy. To compensate the thermal expansion, the pre-stretch of the ball screw is used as the compensation method. At the time of manufacturing, the length of the threaded portion of the ball screw is less than the nominal length (lead by the number of turns) a pre-stretch. The pre-stretch capacity is slightly larger than the thermal expansion. In the assembly process, the lead screw is stretched by a certain tensile structure so that the threaded portion of the lead screw reaches the nominal length. In the service process, the thermal expansion counteracts a part of the pre-stretch capacity, and then the tensile stress of the lead screw decreases, but the length of the screw shaft does not change. The above method ensures that the pitch accuracy is not affected by the thermal expansion. Under no-load conditions, the preload F t applied by the certain tensile structure is equal to the axial load F a . So F a ¼ F t . The greater the preload F t , the more is the heat generation of the screw shaft. Moreover, the greater the pre-stretch capacity, the greater is the preload F t . So the heat generation of the screw shaft increases with the pre-stretch capacity. The change in the pre-stretch capacity leads to the change in the heat generation of the screw shaft and the pre-stretch capacity and the heat generation of the screw shaft are positively correlated. The ball screw is fixed at one end and an axial force is applied to the other end. When the temperature rises, the lead screw generates the compressive stress. The thermal stress distributes in the longitudinal direction, and the induced elongation is expressed as
DL2 ¼ Lr=E
ð64Þ
where r denotes the stress and r ¼ P=A, wherein A and P denote the cross-section area and the axial force, respectively. So the actual elongation is expressed as
Z
L
DL1 ¼ DLðx; tÞ þ DL2 ¼ DLðx; t Þ ¼
ða½T ðx; tÞ T 0 Þdx Lr=E
ð65Þ
0
The pre-stretch capacity is expressed as
DL3 ¼ a D T L
ð66Þ
where D T denotes the average temperature rise. Under no-load conditions, the preload F t applied by the preload nut is equal to the axial load F a . So F a ¼ F t . The positioning error of the feed drive system is expressed as Eq. (67). To ensure that the ball screw feed drive system works properly, it
33
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538 Table 15 Temperature rises and thermal errors. Temperature rise(°C)
Thermal errors(mm)
Pre-stretch capacities (°)
Moving nut
Front bearing housing
Rear rearing housing
Maximum
0 120 150 180
3.92 3.88 4.04 3.80
4.48 4.64 5.11 5.86
3.79 3.91 4.07 5.00
43.7 41.1 32.3 30.2
is necessary to ensure that the difference between the actual elongation of the screw shaft and the pre-stretch capacity is not less than zero. Namely,
DL ¼ DL1 DL3 P 0
ð67Þ
According to Eq. (67), the increase in the pre-stretch capacity leads to the decrease in the positioning error of the linear axis with the fixed-support installation. The pre-stretch capacities are 0°, 120°, 150°, and 180° and then the thermal characteristic tests were carried out. Table 15 shows the temperature rises of the critical points and the maximum thermal errors under various pre-stretch capacities. It is seen that the larger the pre-stretch capacity, the smaller is the maximum thermal error of the feed drive system, and the better is the precision stability of the machine tool. According to Eq. (67), the positioning error of the feed drive system will be reduced with the increase in the pre-stretch capacity, which is consistent with that shown in Table 15. However, it should be noted that the pre-stretch capacities should not be as large as possible, and the excessive pre-stretch capacity will cause the excessive temperature rise of bearings. If the temperature rise is too high, the wear will intensify and even seizure happens to the bearing. Moreover, it is seen that the larger the pre-stretch capacity, the larger the temperature rise of the bearing housing because the increase in the pre-stretch capacity leads to an increase in the bearing load, and the increase in the friction heat causes the increase in the temperature of the bearing housing. In different pre-stretch states, the temperature rise of the nut does not change much because the pre-stretch capacity of the lead screw is still relatively small, and then the effect of the pre-stretch capacity on the fit between the nut and the screw is small. It is revealed that the effect of the pre-stretch of the lead screw on the thermal error of linear axes is realized by changing the deformation of the screw shaft. The pre-stretch can compensate for part of the thermal elongation of the screw shaft. The effect of the pre-stretch of the lead screw on the temperatures and thermal error of linear axes is significant. Moreover, the greater the pre-stretch of the lead screw, the higher is the temperature of linear axes, and the smaller the thermal error of linear axes. To balance the temperature and thermal error, a reasonable pre-stretch capacity should be chosen. 6. Conclusions The data-driven thermally-induced error compensation method was proposed for the five-axis machine tools based on the homogeneous transformation method. The compensation component was obtained by analyzing the error transmission chain of five-axis machine tools and was expressed as the matrix of the error terms of the linear axes and spindle system. From the view of error generation mechanism, the thermally-induced error of the linear axis was formulated as the product of the polynomial function with time as its independent variable and the polynomial function with position as its independent variable. Then the thermal and geometric errors were decomposed and identified from the measured positioning error. The thermal error of the spindle system was expressed as the MLRA model of the critical temperatures. Finally, the error compensation was carried to verify the effectiveness of the data-driven thermally-induced error compensation method. The following conclusions are obtained from the previous analysis. (1) The systematic data-driven thermally-induced error compensation method of the five-axis machine tools is proposed based on the homogeneous transformation method. Then the relationship between the compensation components and the thermal error components in the linear axes and spindle system is revealed by the differential movement of the compensation axis. The thermally-induced error model is robust and effective at different thermal effects, and the data-driven thermally-induced error compensation method could reduce by the machining error more than 85% and 37% with the present error compensation compared with that without compensation and that without traditional error compensation, respectively. (2) The relationship between the thermal expansion, the thermal error, and the resultant positioning error of the linear axis is explored. The temperature of the screw shaft is regarded as the intermediate variables to simplify the decomposition of the thermally-induced error from the measured positioning error of linear axes. For the linear axis with fixed-support installation, the thermal error is formulated as the product of the polynomial function with time as its independent variable and the polynomial function with position as its independent variable. With the aid of the
34
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
polynomial function, the thermal error modeling of the linear axis becomes ever simple and easy. The thermal error is the main component of the positioning error of the linear axis. The decomposition of the thermal error from the positioning error is realized by the polynomial fitting. (3) The direct and accurate calculation of the thermal error of the linear axis is feasible from the temperature rise. The temperature response of the screw shaft is the superposition of three responses of equivalent temperatures. The thermal error is the integral of the temperature rise on the length of the lead screw. It is revealed that the MLRA method is effective to establish the correlation between the thermally-induced error and critical temperatures of the spindle system. The dependence of the thermal error on the feed rate and rotational speed is strong. Both the temperature rise and thermally-induced error of the spindle system increase with the rotational speed and feed rate. The lag effect of the thermally-induced error is obvious respect to the fast change in the feed rate and rotational speed. (4) The effect of the heat source on the thermal error is realized by changing the temperatures of the heat sources. Then the heat entering linear axes is changed. By substituting the change in the heat generation into the thermal field and thermal error models of linear axes, the temperature and thermal error are predicted accurately. The effect of the prestretch of the lead screw on the thermal error of linear axes is realized by changing the deformation of the screw shaft. The pre-stretch can compensate for part of the thermal elongation of the screw shaft. The effect of the pre-stretch of the lead screw on the temperatures and thermal error of linear axes is significant. Moreover, the greater the prestretch of the lead screw, the higher is the temperature of linear axes, and the smaller the thermal error of linear axes. To balance the temperature and thermal error, a reasonable pre-stretch capacity should be chosen.
Acknowledgments This research was supported by the National Natural Science Foundation of China (51905057), the Natural Science Foundation Project of Chongqing, Chongqing Science and Technology Commission (cstc2019jcyj-msxmX0050), and the Fundamental Research Funds for the Central Universities (2018CDXYJX0019). Declaration of Competing Interest We declare that we have no financial and personal relationships with other people or organizations that can inappropriately influence our work; there is no professional or other personal interest of any nature or kind in any product or company that could be construed as influencing the position presented in, or the review of, the manuscript. References [1] J. Mayr, J. Jedrzejewski, E. Uhlmann, et al, Thermal issues in machine tools, CIRP Ann.-Manuf. Technol. 61 (2) (2012) 771–791. [2] Y. Altintas, M.R. Khoshdarregi, Contour error control of CNC machine tools with vibration avoidance, CIRP Ann.-Manuf. Technol. 61 (1) (2012) 335–338. [3] C. Ma, X. Mei, J. Yang, et al, Thermal characteristics analysis and experimental study on the high-speed spindle system, Int. J. Adv. Manuf. Technol. 79 (1–4) (2015) 469–489. [4] C. Ma, L. Zhao, X. Mei, et al, Experimental and simulation study on the thermal characteristics of the high-speed spindle system, ARCHIVE Proc. Inst. Mech. Eng. Part C J. Mech. Eng. Sci. 1989-1996 203–210 (6) (2016) 6. [5] C. Ma, J. Yang, L. Zhao, et al, Simulation and experimental study on the thermally induced deformations of high-speed spindle system, Appl. Therm. Eng. 86 (2015) 251–268. [6] J. Liu, C. Ma, S. Wang, et al, Thermal-structure interaction characteristics of a high-speed spindle-bearing system, Int. J. Mach. Tools Manuf. (2019), https://doi.org/10.1016/j.ijmachtools.2018.10.004. [7] J. Liu, C. Ma, S. Wang, et al, Thermal contact resistance between bearing inner ring and shaft journal, Int. J. Therm. Sci. 138 (2019) 521–535. [8] J. Liu, C. Ma, S. Wang, Thermal contact conductance between rollers and bearing rings, Int. J. Therm. Sci. (2020,), https://doi.org/10.1016/j. ijthermalsci.2019.106140. [9] C. Ma, L. Zhao, H. Shi, et al, A geometrical-mechanical-thermal predictive model for thermal contact conductance in vacuum environment, Proc. Inst. Mech. Eng., Part B: J. Eng. Manuf. 230 (8) (2016) 1451–1464. [10] A.L. Sweet, J.F. Tu, Tolerance design for the fit between bore and shaft for precision assemblies with significant error-scaling problems, Int. J. Prod. Res. 45 (22) (2007) 5223–5241. [11] C. Xia, J. Fu, J. Lai, et al, Conjugate heat transfer in fractal tree-like channels network heat sink for high-speed motorized spindle cooling, Appl. Therm. Eng. 90 (2015) 1032–1042. [12] J. Liu, C. Ma, S. Wang, et al, Contact stiffness of spindle-tool holder based on fractal theory and multi-scale contact mechanics model, Mech. Syst. Sig. Process. 119 (2019) 363–379. [13] J. Liu, C. Ma, S. Wang, Precision loss modeling method of ball screw pair, Mech. Syst. Sig. Process. (2020,), https://doi.org/10.1016/j.ymssp.2019.106397. [14] E. Creighton, A. Honegger, A. Tulsian, et al, Analysis of thermal errors in a high-speed micro-milling spindle, Int. J. Mach. Tools Manuf. 50 (4) (2010) 386–393. [15] H. Zhao, J. Yang, J. Shen, Simulation of thermal behavior of a CNC machine tool spindle, Int. J. Mach. Tools Manuf. 47 (6) (2007) 1003–1010. [16] J.J. Kim, Y.H. Jeong, D.W. Cho, Thermal behavior of a machine tool equipped with linear motors, Int. J. Mach. Tools Manuf. 44 (7) (2004) 749–758. [17] B. Bossmanns, J.F. Tu, A power flow model for high speed motorized spindles—heat generation characterization, J. Manuf. Sci. Eng. 123 (3) (2001) 494– 505. [18] B. Bossmanns, J.F. Tu, A thermal model for high speed motorized spindles, Int. J. Mach. Tools Manuf. 39 (9) (1999) 1345–1366. [19] C.H. Chien, J.Y. Jang, 3-D numerical and experimental analysis of a built-in motorized high-speed spindle with helical water cooling channel, Appl. Therm. Eng. 28 (17) (2008) 2327–2336. [20] J. Liu, C. Ma, S. Wang, et al, Thermal boundary condition optimization of ball screw feed drive system based on response surface analysis, Mech. Syst. Sig. Process. 121 (2019) 471–496. [21] Y. Kang, C.W. Chang, Y. Huang, et al, Modification of a neural network utilizing hybrid filters for the compensation of thermal deformation in machine tools, Int. J. Mach. Tools Manuf. 47 (2) (2007) 376–387.
J. Liu et al. / Mechanical Systems and Signal Processing 138 (2020) 106538
35
[22] J.P. Kruth, P. Vanherck, V.D.B. Christophe, Compensation of static and transient thermal errors on CMMs, Ann. CIRP 50 (1) (2001) 377–380. [23] H. Yang, J. Ni, Dynamic neural network modeling for nonlinear, nonstationary machine tool thermally induced error, Int. J. Mach. Tools Manuf. 45 (4) (2005) 455–465. [24] R. Ramesh, M.A. Mannan, A.N. Poo, Error compensation in machine tools — a review: Part II: thermal errors, Int. J. Mach. Tools Manuf. 40 (9) (2000) 1257–1284. [25] S.H. Yang, K.H. Kim, K.P. Yong, Measurement of spindle thermal errors in machine tool using hemispherical ball bar test, Int. J. Mach. Tools Manuf. 44 (2–3) (2004) 333–340. [26] W.T. Lei, Y.Y. Hsu, Accuracy enhancement of five-axis CNC machines through real-time error compensation, Int. J. Mach. Tools Manuf. 43 (9) (2003) 871–877. [27] Y. Zhang, J. Yang, H. Jiang, Machine tool thermal error modeling and prediction by grey neural network, Int. J. Adv. Manuf. Technol. 59 (9–12) (2012) 1065–1072. [28] J. Yang, J. Yuan, J. Ni, Thermal error mode analysis and robust modeling for error compensation on a CNC turning center, Int. J. Mach. Tools Manuf. 39 (9) (1999) 1367–1381. [29] Y. Hong, J. Ni, Adaptive model estimation of machine-tool thermal errors based on recursive dynamic modeling strategy, Int. J. Mach. Tools Manuf. 45 (1) (2005) 1–11. [30] C. Ma, L. Zhao, X. Mei, et al, Thermal error compensation based on genetic algorithm and artificial neural network of the shaft in the high-speed spindle system, Proc. Inst. Mech. Eng., Part B: J. Eng. Manuf. 231 (5) (2017) 753–767. [31] C. Ma, L. Zhao, X. Mei, et al, Thermal error compensation of high-speed spindle system based on a modified BP neural network, Int. J. Adv. Manuf. Technol. 89 (9–12) (2017) 3071–3085. [32] Z.Z. Xu, C. Chang, L.J. Liang, et al, Study on a novel thermal error compensation system for high-precision ball screw feed drive (2nd report: Experimental verification), Int. J. Precis. Eng. Manuf. 16 (10) (2015) 2139–2145. [33] W.S. Yun, S.K. Kim, W.C. Dong, Thermal error analysis for a CNC lathe feed drive system, Int. J. Mach. Tools Manuf. 39 (7) (1999) 1087–1101. [34] H. Shi, B. He, Y. Yue, et al, Cooling effect and temperature control of oil cooling system for ball screw feed drive system of precision machine tool, Appl Therm Eng. 161 (10) (2019) 114150. [35] Z.Z. Xu, X.J. Liu, S.K. Lyu, Study on positioning accuracy of nut/shaft air cooling ball screw for high-precision feed drive, Int. J. Precis. Eng. Manuf. 15 (1) (2014) 111–116. [36] H. Pahk, S.W. Lee, Thermal error measurement and real time compensation system for the CNC machine tools incorporating the spindle thermal error and the feed axis thermal error, Int. J. Adv. Manuf. Technol. 20 (7) (2002) 487–494. [37] A.C. Okafor, ERTEKIN, et al, Derivation of machine tool error models and error compensation procedure for three axes vertical machining center using rigid body kinematics, Int. J. Mach. Tools Manuf. 40 (8) (2000) 1199–1213. [38] Z. Li, J. Yang, K. Fan, et al, Integrated geometric and thermal error modeling and compensation for vertical machining centers, Int. J. Adv. Manuf. Technol. 76 (5–8) (2015) 1139–1150. [39] Y. Wang, G. Zhang, K.S. Moon, et al, Compensation for the thermal error of a multi-axis machining center, J. Mater. Process. Technol. 75 (1–3) (1998) 45–53. [40] J.H. Jung, J.P. Choi, S.J. Lee, Machining accuracy enhancement by compensating for volumetric errors of a machine tool and on-machine measurement, J. Mater. Process. Technol. 174 (1–3) (2006) 56–66. [41] Y. Lu, A new approach to thermally induced volumetric error compensation, Int. J. Adv. Manuf. Technol. 62 (9–12) (2012) 1071–1085. [42] Incropera, Frank P. Fundamentals of Heat and Mass Transfer. Wiley, 1985:139–162. [43] ISO 230-3, Test Code for Machine Tools Part 3: Determination of Thermal Effects, ISO copyright office, Switzerland, 2007.. [44] H. Shi, C. Ma, J. Yang, et al, Investigation into effect of thermal expansion on thermally induced error of ball screw feed drive system of precision machine tools, Int. J. Mach. Tools Manuf. 97 (2015) 60–71.