Int J Adv Manuf Technol DOI 10.1007/s00170-016-9924-2

ORIGINAL ARTICLE

An optimal convective heat transfer coefficient calculation method in thermal analysis of spindle system Feng Tan 1 & Qin Yin 1 & Guanhua Dong 1 & Luofeng Xie 1 & Guofu Yin 1

Received: 19 July 2016 / Accepted: 18 December 2016 # Springer-Verlag London 2017

Abstract Accurate calculation of convective heat transfer coefficients (CHTCs) is the premise to obtain accurate results of finite element thermal analysis (FETA). This paper presents an optimization method of CHTCs in steady state FETA of a spindle system using genetic algorithm (GA), which is especially suitable for complicated FETA with many CHTCs, such as 10 CHTCs. Firstly, the temperature values at thermal critical points, the temperature field, and the thermal fluctuations of a spindle system were obtained by experiment. Secondly, based on similarity theory, the initial values of CHTCs, which are inaccurate but close to real ones, were calculated by empirical formulas, and then, the FETA of the spindle system was conducted. The inaccurate initial CHTCs lead to inaccurate FETA results, so the optimal CHTCs closer to real ones were searched by GA automatically. The CHTCs were regarded as the target values of interest, and the root-mean-square error (RMSE) between the experimental temperature values and the simulation temperature values was regarded as the fitness function. After 200 generations, the optimal CHTCs were searched with which the RMSE was the lowest. Finally, by conducting the FETA of the spindle system with the searched optimal CHTCs and comparing with the experimental results, the effectiveness of the proposed method was proved.

* Guofu Yin [email protected] 1

School of Manufacturing Science and Engineering, Sichuan University, Chengdu 610065, People’s Republic of China

Keywords Convective heat transfer coefficient . Finite element thermal analysis . Genetic algorithm . Spindle system

1 Introduction With the development of modern manufacturing technology, there is a growing demand for higher machining precision machine tools. Numerous studies have shown that thermal errors are becoming a major factor affecting the machining accuracy of machine tools and account for as high as 70% of total errors affecting the accuracy of machine tools [1–3]. Besides, the spindle system (motorized spindle and spindle box) is the biggest heat source of the machine tool, so are the major source of machine tool thermal errors [4, 5]. Therefore, it is necessary to analyze and optimize the thermal characteristics of the machine tool to reduce the thermal error and improve machining precision, especially the spindle system. Nevertheless, the key of analysis and optimization is to establish an accurate thermal characteristics model of the machine tool. Usually, experimental method and finite element (FE) method are the two main methods for thermal characteristics modeling. Although the thermal characteristics of the spindle system can be obtained through experimental method accurately, this method has some disadvantages such as being time consuming and costly to conduct the experiment and hard to acquire comprehensive thermal characteristics of the spindle system. More importantly, the experiment can only be conducted after a machine tool is manufactured. On the contrary, the FE method can overcome the shortcomings above-mentioned and can analyze the comprehensive thermal characteristics of a spindle system so as to conduct the optimal design of the next generation spindle system at its design stage. Nevertheless, in most cases, FE simulation results differ greatly from the actual values because the

Int J Adv Manuf Technol

boundary conditions in finite element thermal analysis (FETA) are often not consistent with actual situation. Therefore, it is essential to calculate the accurate boundary conditions in FETA. The accuracy of FETA is most likely affected by thermal sources, thermal contact resistances (TCRs), and convective heat transfer coefficients (CHTCs) [6, 7]. For the spindle system, the major thermal sources are the front and rear bearings and the built-in motor, and the heat power can be calculated by mature methods [8, 9]. The TCRs between different parts can be obtained experimentally or theoretically [10–13]. The CHTCs reflect the capability of heat taken away by the surrounding fluid from the spindle system surfaces, and according to fluid flow state, the convective heat transfer can be divided into forced convective heat transfer and natural convective heat transfer [14]. Moreover, the CHTCs are influenced by many factors such as fluid properties, shape of the heat transfer surfaces, fluid velocity, and temperature differences between the heat transfer surfaces and fluid; thus, the CHTCs are usually calculated by similarity theory but with low accuracy. Many researchers have studied the method to calculate more accurate CHTCs to improve the accuracy of FETA. Min et al. [15] simulated the temperature field of a ball screw system but with low accuracy initially. So, they manually changed the CHTC from the initial value with a step value to a certain value, with which the simulated temperature agreed well with the experimental temperature. Yang et al. [16] constructed a neural network with temperatures as input and CHTCs as output to iteratively search the optimal CHTCs in FETA of a ball screw system. Chow et al. [17] compared the FETA results of a carriage with those obtained through experiment and pointed out that it is important to know the thermal boundary conditions including CHTCs to accurately conduct the FETA. Li et al. [18] established the response surface model between CHTCs and simulation error and searched the best CHTCs with the lowest simulation error. According to the above literature, the current optimization method of CHTCs are unsatisfactory especially not suitable for complicated FETA with many CHTCs, thus resulting in poor simulation accuracy. To solve the problem, we proposed a method using GA to calculate the optimal CHTCs automatically and iteratively in complicated FETA. There are many convective heat transfers occurring during the rotation of the spindle system and they vary greatly. So, we take a steady state FETA of a spindle system of a horizontal machining center as an example. Section 2 described the thermal experiment of the spindle system. Then, Sect. 3 conducted the FETA of the spindle system with the initially calculated CHTCs using empirical formulas. Because of the big difference between experimental results and simulation results,

Sect. 4 demonstrated the optimization method of the CHTCs of the spindle system using GA. Afterwards, the effectiveness of the proposed method was proved in Sect. 5. Finally, conclusions were drawn in Sect. 6.

2 Thermal experiment of the spindle system The spindle system of a horizontal machining center THM6380 is shown in Fig. 1, consisting of a spindle box and a motorized spindle. As detailed in Fig. 2, the spindle system mainly consists of a spindle box, a shaft, four front bearings (angular contact ball bearing), one rear bearing (cylindrical roller bearing), spacers, bearing sleeves, spindle sleeve, water jacket, end caps, and stator and rotor of built-in motor. The thermal experiment setup of the spindle system is depicted in Fig. 3. To obtain an overall temperature field distribution of the spindle system, seven platinum resistance temperature sensors (T1–T7) were fixed on the seven thermal critical points along the axial direction of the spindle system as the temperature gradient varies almost along the axial direction. T1–T4 were used to measure the temperatures of the front side of the motorized spindle, T5–T7 were used to measure the temperatures of the rear side of the motorized spindle, and the environmental temperature was measured by T8. The real-time temperatures were recorded every 5 min. Meanwhile, a thermal imager was used to capture the temperature field of the front part of the spindle system also at 5-min intervals, which could be used to compare with the simulated temperature field visually. To measure the deformations of the spindle system, a test bar was assembled in the spindle system. As shown in Fig. 3, three capacitive displacement transducers were fixed on the worktable to measure the thermal deformations of the test bar in the X, Y, and Z directions, and the deformations were also recorded every 5 min in real time. The spindle speed was 3000 rpm, and in order to obtain the real thermal characteristics of the spindle system, the water cooling system was turned off. Theoretically, it will take infinite time to reach the thermal steady state. However, in practice, when the currently recorded temperatures were almost unchanged with respect to those recorded at last moment, we take the current state as the steady state of the spindle system [19]. And it took about 5 h to reach that relatively steady state. The experimental results show that the steady state temperatures of T1–T7 were 30.9, 38.7, 39.3, 38.8, 42.8, 39.6, and 40.5 °C, respectively. The temperature of T1 was the lowest occurred on the front-end cap with the reason that a large amount of heat was taken away by the rotating test bar and the big TCR between the front-end cap with its neighbor parts. The temperature of T5 was the highest occurred on the rear side of the motorized spindle because of close to the rear bearing and the poor heat dissipation condition on the rear

Int J Adv Manuf Technol Fig. 1 The 3D model of the spindle system

side. T2–T4 have relative high temperatures because they are close to the front bearings and the built-in motor. T5–T7 have relative high temperatures because they are close to the builtin motor and the rear bearing. The temperatures of T2 and T4 and the temperatures of T3 and T6 were virtually identical in value, respectively. The ambient temperature was stable at about 26.7 °C, and it was regarded as the initial temperature value in FETA. In addition, the thermal deformations in X and Y directions of the spindle system were almost zero, and the thermal elongation in Z direction was 137.1 μm. The steady state thermal image of the front side of the spindle system is shown in Fig. 4. As shown, the steady state temperature field distribution of the front part of the spindle system is almost uniform. The three square regions marked with A, B, and C in Fig. 4 were the emissivity-adjusted areas because of the high specular surfaces.

Fig. 2 Structure of the spindle system

3 Steady state FETA of the spindle system The steady state FETA process included three steps: (1) establishing the FE model of the spindle system; (2) calculating the boundary conditions, including heat sources, TCRs, and CHTCs; and (3) executing the FETA. 3.1 Establishing the FE model The 3D model of the spindle system is shown in Fig. 1, comprising a motorized spindle and a spindle box. Small details not impacting on analysis such as holes, fillets, and chamfers were ignored. The spindle system was made from several main materials, and these materials’ properties were assigned to their corresponding parts in ANSYS Workbench simulation software. The material assignment details and the material

Int J Adv Manuf Technol

two major heat sources: (1) friction heat generated by front and rear bearings and (2) heat generated by the built-in motor. (1) Friction heat generated by front and rear bearings The heat generated by a bearing can be calculated [8] by H b ¼ 1:047  10−4 nM

ð1Þ

where Hb is the heat power (W), n is the bearing rotation speed (rpm), and M is the total frictional torque (N mm) of the bearing. The total frictional torque mainly consists of two parts: M ¼ Mυ þ Ml

ð2Þ

where Mυ is the viscous friction torque (N mm) caused by the lubricant viscosity and Ml is the mechanical friction torque (N mm) caused by the applied load. The calculation of Mυ is related to the value of υ0n, and in this study υ0n ≥ 2000, so the viscous friction torque Mυ can be calculated by Fig. 3 Thermal experiment setup of the spindle system

properties are shown in Tables 1 and 2, respectively. The motorized spindle was meshed into hexahedral mesh, and the spindle box was meshed into tetrahedral mesh with a total of 317,340 nodes and 276,436 elements, shown in Fig. 5.

M υ ¼ 10−7 f 0 ðυ0 nÞ2=3 d 3m

ð3Þ

where f0 is a coefficient related to the bearing type and lubricant method, υ0 is the kinematic viscosity (mm2/s) of the lubricant, and dm is the mean diameter (mm) of the bearing. The mechanical friction torque Ml can be calculated by

3.2 Calculating the boundary conditions

M l ¼ f 1 F β dm

3.2.1 Calculating the heat sources

where f1 is a coefficient depending on the bearing design and relative bearing load and Fβ is the equivalent bearing load (N). The detailed calculation formulas of f1 and Fβ can be referred to reference [8].

The idling thermal characteristics of the spindle system were simulated in this study. During idling, the spindle system has Fig. 4 Steady state thermal image of the front side of the spindle system

ð4Þ

Int J Adv Manuf Technol Material assignment details

Table 1 Material

Corresponding parts

45 GCr15

Spacers, nuts, end caps, sleeves, water jacket, stator, rotor, spindle nose Bearings

20CrMnMoH HT250

Shaft Spindle box

used in this study were obtained through an experimental method [11]. In actual FETA of the spindle system, the thermal contact conductance (TCC), which is the reciprocal of the TCR, was assigned to the contact surfaces, so the TCC values of key contact surfaces were calculated and are shown in Table 3.

3.2.3 Calculating the CHTCs Through computation, the heat power of each front and rear bearing are 65.4 and 39.6 W, respectively. (2) Heat generated by built-in motor The copper loss, iron loss, and other additional losses of built-in motor are all converted to heat power. It is hard to accurately calculate all the losses, so usually the efficiency analysis method was used to calculate the heat power of the built-in motor [9]. The approximate calculation equation of the heat power generated by stator and rotor of the built-in motor can be expressed as Hm ¼

During the air cutting of the spindle system, there exists four kinds of convective heat transfer: (1) forced convective heat transfer between rotating surfaces and the surrounding air, (2) forced convective heat transfer between rotor and stator of the built-in motor, (3) forced convective heat transfer between oil mist and bearings, and (4) natural convective heat transfer between stationary surfaces and the ambient air. The calculation of CHTCs is a complex function influenced by many factors, as follows: h ¼ f ðv; D; λ; μ; c; ρ; x; y; z; φ⋯Þ

1−ηm nM ηm 9550

ð5Þ

where ηm is the mechanical efficiency of the built-in motor and ηm = 0.9 in this study; M is the total frictional torque caused by the front and rear bearings as calculated by Eq. (2). The calculated heat power of the built-in motor is 66.9 W. 3.2.2 Calculating the TCRs In actual situation, because of the existence of the TCRs between contact components, when heat flows from a hotter component to a colder component, there is an obvious temperature drop between two contact surfaces. Factors affecting TCRs are very complex, including the materials of contact components, the filler between the contact surfaces, surface roughness, and the contact pressure [10]. The TCRs between different parts greatly affect the heat transfer. In this study, TCRs at key contact surfaces of the spindle system were considered: (1) between the bearing inner ring and the shaft, (2) between the bearing outer ring and the spindle sleeve, and (3) between the front-end cap and its neighbor parts. The TCRs

ð6Þ

where h is the CHTC; v is the fluid flow velocity; D is the feature size; λ is the thermal conductivity of the fluid; μ is the dynamic viscosity of the fluid; c is the specific capacity of the fluid; ρ is the density of the fluid; x, y, and z are the space coordinates; and φ is the heat transfer surface geometry. Because of the complexity of the calculation, these CHTCs are generally calculated by empirical equations based on similarity theory, thus resulting in an inaccurate simulation results. In actual situation, the convective heat transfers between the spindle system and the surrounding fluids are complicated. In order to improve the FETA accuracy, 10 CHTCs were defined according to the convective heat transfer state and the geometry structure of the spindle system, as shown in Fig. 6. In order to differentiate, forced CHTCs are marked with red and orange lines alternately and natural CHTCs are marked with blue and green lines alternately. (1) Forced CHTCs between rotating surfaces and surrounding air When the shaft is rotating, relative movement occurs between the rotating surfaces and the air, thus producing forced

Table 2 Material properties Materials

Young’s modulus (GPa)

Poisson’s ratio

Density (kg/m3)

Coefficient of thermal expansion (10−5 m/K)

Thermal conductivity (W/(m K))

Specific heat capacity (J/(kg K))

45 GCr15 20CrMnMoH HT250 Air

209 219 207 138 –

0.269 0.3 0.254 0.156 –

7890 7830 7870 7280 1.293

1.17 1.20 1.12 0.82 –

48 44 44 45 0.02454

450 460 460 510 1.005

Int J Adv Manuf Technol

where Re is the Reynolds number and Pr is the Prandtl number of the flowing air, which can be calculated by the hydrodynamic formulas of the flowing air. According to Eqs. (7)–(9), the forced CHTC between the shaft end surfaces and surrounding air is h1 = 91.6W/(m2 K) and the forced CHTC between the shaft outer surfaces and surrounding air is h2 = 66.2W/(m2 K). (2) Forced CHTC between rotor and stator of the built-in motor When the built-in motor is rotating, forced convective heat transfer occurs between the rotor and the stator. Similar to the forced CHTC between rotating surfaces and surrounding air, the CHTC between the rotor and the stator of the built-in motor can be calculated by [21] h¼

Nuλair Δ  10−3

ð10Þ

where Δ is the gap between the rotor and the stator and Nu is the Nusselt number computed by Fig. 5 Mesh model of the spindle system

Nu ¼ 0:23ðΔ=rÞ0:25 Re0:5

convective heat transfer. There are two kinds of rotating surfaces, including end surfaces and outer surfaces. For the end surfaces, the forced convective heat transfer coefficient can be calculated by  pffiffiffiffiffiffiffiffiffiffiffi ð7Þ h ¼ 28 1 þ 0:45u where u is the ambient air velocity (m/s) around the end surfaces, which can be calculated by a simple linear velocity formula. For the outer surfaces, the forced convective heat transfer coefficient can be calculated by the following Nusselt’s equation [20]: h¼

Nuλair d o  10−3

ð8Þ

where λair is the thermal conductivity of the air (W/(m K)); do is the diameter (mm) of the outer surfaces; and Nu is the Nusselt number calculated by Nu ¼ 0:133Re2=3 Pr1=3

Table 3

ð9Þ

ð11Þ

where r is outer radius of the rotor (mm) and Re is the Reynolds number computed by the hydrodynamic formulas of the flowing air between the rotor and the stator. Based on Eqs. (10) and (11), the CHTC between the rotor and the stator is h3 = 65.3W/(m2 K). (3) Forced CHTC between oil mist and bearings The motorized spindle bearings are lubricated by oil mist ejected from the nozzle. Apart from lubrication, the oil mist lubricant can also take away a certain amount of heat. The forced CHTC between oil mist lubricant and bearings were determined by experiment [22], and the typical value 42 W/ (m2 K) was adopted in this study since the rotation speed is 3000 r/min. That is to say, h4 = h5 = 42 W/(m2 K). (4) Natural CHTC between stationary surfaces and the ambient air All other stationary surfaces exposed to ambient air, such as the spindle box and the rear part of the motorized spindle, dissipate heat through natural convective heat transfer. According to previous experience, 9.7 W/(m2 K) was taken as the natural CHTC [22]. That is to say, h6 = h7 = h8 = h9 = h10 = 9.7 W/(m2 K). The calculated initial values of all CHTCs through empirical formulas are h0 = [91.6, 66.2, 65.3, 42, 42, 9.7, 9.7, 9.7, 9.7, 9.7].

TCC values of key contact surfaces (W/(m2 K))

Contact surfaces

TCC values

Bearing inner ring and shaft Bearing outer ring and spindle sleeve Front-end cap and its neighbor parts

6850 1490 66

3.3 Thermal simulation of the spindle system In order to simulate the steady state temperature field of the spindle system, the steady state FETA was performed. The thermal experiment of the spindle system was conducted in a relatively constant temperature workshop, so during the

Int J Adv Manuf Technol Fig. 6 The defined 10 CHTCs of the spindle system

steady sate FETA of the spindle system, the boundary conditions, i.e., heat sources, TCCs, and CHTCs, were assumed to be constant. The initial temperature value was set to be 26.7 °C, which is identical to the ambient temperature. The material properties and the calculated boundary conditions were applied to the FE model. Figure 7a shows the simulated temperature field of the spindle system and temperature values of key temperature points T1–T7. Table 4 shows the comparison of the experimental and simulated temperature values of the key temperature points (T1–T7), and the maximum simulation error is about 46.3%. It is obviously seen that there exist distinct differences between the experimental and simulated temperatures. By applying the simulated temperature field as body load, and imposing fixed constraints on the surfaces the sliders installed, the thermal fluctuations of the spindle system were simulated. Figure 7b shows the simulated axial (Z direction) thermal fluctuations of the spindle, and the value of the spindle nose is about 151.28 μm with a simulation error about 10.34%. Comparing Fig. 7a with Fig. 4, the simulated and experimental temperature fields of the front side of the spindle system are also very different. So, this study presents a method that can be used to improve the steady state FETA accuracy of the spindle system from the perspective of optimizing the CHTCs.

target values of interest, i.e., the chromosome in GA, and the RMSE between the simulation temperature values and the experimental temperature values is regarded as the fitness

(a)

4 Searching the optimal CHTCs using genetic algorithm The GA can search the global minimum with the merits of self-organization, self-adaption, and self-learning. So, the CHTCs are optimized by GA, which is widely used in variable optimization [23–27]. The 10 CHTCs are regarded as the

(b) Fig. 7 a Simulated temperature field and temperature values of key temperature points T1–T7. b Simulated axial (Z direction) thermal fluctuations

Int J Adv Manuf Technol Table 4

Comparison of experimental and simulated temperature values

Temperature points

Experimental/ °C

Simulated/ °C

Error/ °C

Error/ %

T1

30.9

45.17

14.27

46.3

T2

38.7

50.12

11.42

29.5

T3 T4

39.3 38.8

44.76 36.87

5.46 −1.93

14 −4.9

T5

42.8

33.24

−9.56

−22.4

T6 T7

39.6 40.5

32.48 31.73

−7.12 −8.77

−17.9 −21.7

function. When the fitness value is less than a certain value or the generation reaches a certain number, the GA stops, and the target values of interest of the current generation are regarded as the optimal CHTCs. The optimization flowchart is shown in Fig. 8, mainly including the following operations: 4.1 Coding operation The target values of interest length are set as the number of CHTCs of 10. And based on the initial CHTCs calculated by empirical formulas, the optimization bounds can be determined. The lower and upper bounds are one tenth and six times of the initial CHTCs, respectively, to ensure a large searching scope, as shown in Table 5. Since the CHTCs are directly used in the FETA, the real numbers of the CHTCs are adopted as the coding contents. An initial population with a population size of 10 is generated by coding the 10 CHTCs within the optimization bounds.

(a)

(b)

4.2 Fitness function The purpose of this study was to improve the steady state FETA accuracy of the spindle system, that is to say reducing the differences between the simulation temperatures and experimental temperatures. First, the corresponding 10 CHTCs, i.e., the target values of interest, were assigned to the FE model. Then, the FETA with the assigned CHTCs was conducted in ANSYS Workbench simulation software. Finally, the RMSE between the experimental temperature values and the simulation temperature values was regarded as the fitness value, which can be expressed as follows: sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 7  s e 2 F¼ ð23Þ ∑ T −T 7 i¼1 i i

Fig. 8 a Optimization flowchart using GA. b Data communication between MATLAB and ANSYS

where T si (i = 1 , 2 , ⋯ , 7) are the simulation temperature values and T ei (i = 1 , 2 , ⋯ , 7) are the experimental temperature values obtained in thermal experiment in Sect. 2.

pi ¼

4.3 Selection operation There are many selection methods in selection operation, such as roulette selection method and championship selection method. In this study, the roulette selection method was adopted, namely, selecting individuals based on fitness ratio. The probability of each individual can be selected can be expressed as . ð24Þ f i ¼ 1 Fi fi

ð25Þ

N

∑ f j¼1

j

Int J Adv Manuf Technol Table 5 Optimization bounds of CHTC values (W/(m2 K))

CHTCs

h2

h1

Lower bound

9.16

Upper bound

h3 6.62

549.6

397.2

where Fi is the fitness value of the ith individual calculated by Eq. (23), pi is the selection probability of the ith individual, and N is the number of individuals in a population. The reciprocal of the fitness value was calculated due to the purpose of obtaining the smaller fitness value. After the selection operation, the individuals with big fitness values are eliminated, while the small ones are selected to heredity to the next generation. At the same time, the average fitness, the best fitness, and the best target values of interest corresponding to the best fitness of each generation were recorded.

4.4 Cross operation The individuals were coded using real numbers, so the real number cross method was adopted. The cross operation in the jth position between the kth target values of interest ak and the lth target values of interest al can be expressed as  akj ¼ akj ð1−bÞ þ alj b ð26Þ alj ¼ alj ð1−bÞ þ akj b where b is the cross probability which is a random number within [0,1] and b = 0.4 in this study.

4.5 Mutation operation The jth gene in the ith target values of interest aij can be mutated as follows:    aij þ aij −amax f ðg Þ; r > 0:2 ð27Þ aij ¼ aij þ amin −aij f ðg Þ; r ≤ 0:2  g 2 f ðgÞ ¼ r2 1− ð28Þ G

h4 6.53

391.8

h5 4.2

252

4.2 252

h2

h3

h4

h5 h6

h7

h8

h9

0.97 58.2

h8 0.97 58.2

h9 0.97 58.2

h10 0.97 58.2

4.6 Data communication We established the GA optimization program in MATLAB platform, and the data communication between MATLAB and ANSYS is depicted in Fig. 8b. During each data communication, the MATLAB first rewrite the CHTCs to the APDL file and then call ANSYS for batch run with the updated APDL file. Soon after, the newly simulated temperatures Ts were written to a TXT file according to the APDL command. Finally, the MATLAB read the newly simulated temperatures Ts to calculate the RMSE and to determine whether to continue or to stop the loop according to whether the loop end condition is reached.

4.7 Optimization results The population size was set as 20, maximum generation as 200, cross probability as 0.4, and mutation probability as 0.2. Then, the optimal CHTCs can be automatically searched using the GA optimization program. Finally, after 200 generations, the 10 optimal CHTC values were obtained. Comparison of the initially calculated and the optimally searched CHTCs are shown in Table 6. The fitness curve is shown Fig. 9.

Table 6 Comparison of the original calculated and optimal CHTCs (W/(m2 K)) h1

0.97 58.2

h7

probability which is a random number within [0,1]; and r = 0.2 in this study.

where amax and amin are the upper and lower bounds of gene aij, respectively; r2 is a random number; g is the current generation; G is the maximum generation; r is the mutation

CHTCs

h6

h10

Original 91.66 66.2 65.3 42 42 9.7 9.7 9.7 9.7 9.7 Optimized 431.6 211.7 7.3 154.5 26 40.1 49.8 7.8 4.7 1 Fig. 9 Fitness curve

Int J Adv Manuf Technol Fig. 10 a Simulated temperature field and temperature values of key temperature points T1–T7. b Simulated axial (Z direction) thermal fluctuations

(a)

(b) 5 Validation of the optimal CHTCs Table 7

In this section, the 10 optimal CHTCs optimized by GA were applied to the FE model of the spindle system to conduct steady state thermal analysis and static structural analysis. The material properties, heat sources, and TCCs were the same with those in Sect. 3. Using the optimal CHTCs, the simulated temperature field of the spindle system and temperature values of key temperature points T1–T7 were shown in Fig. 10a. Table 7 shows the comparison of the experimental and simulated temperature values using the optimal CHTCs, and the maximum simulation error is about 8.1%. The simulated axial (Z direction)

Comparison of experimental and simulated temperature values

Temperature points

Experimental/ °C

Simulated/ °C

Error/ °C

Error/ %

T1 T2 T3 T4 T5 T6 T7

30.9 38.7 39.3 38.8 42.8 39.6 40.5

33.39 41.48 40.07 38.50 40.40 40.21 40.00

2.49 2.78 0.77 −0.3 −2.4 0.61 −0.5

8.1 7.2 2 −0.8 −5.6 1.5 −1.2

Int J Adv Manuf Technol

thermal fluctuations of the spindle was shown in Fig. 10b, and the value of the spindle nose is about 133.81 μm with a simulation error about −1.2%. Compared with the simulation errors in Sect. 3, the simulation errors in this section using the optimal CHTCs were decreased significantly. Moreover, comparing Fig. 10a with Fig. 4, the simulated and experimental temperature fields of the front side of the spindle system are also very similar. In general, the simulation results using the optimal CHTCs obtained by GA agreed well with the experimental results.

Acknowledgments This research was financially supported by the Key National Science and Technology Projects of China (grant no. 2015BAF27B01) and Science and Technology Support Plan Project of Sichuan Province (grant no. 2014GZ0125).

References 1.

2.

6 Conclusions This paper presented a method using GA for searching the optimal CHTCs in steady state FETA of a spindle system to improve the simulation accuracy. 1. The thermal experiment of the spindle system was conducted. Temperatures at seven thermal critical points were measured using temperature sensors, and the temperature field of the front side of the spindle system was captured by a thermal imager. The thermal fluctuations of the spindle system were measured using capacitive displacement transducers. 2. The steady state FETA of the spindle system was conducted using the initial CHTCs calculated by empirical formulas. However, the simulation results differ greatly from the experimental results with a maximum temperature error up to 46.7%. 3. By regarding the CHTCs as the target values of interest and the RMSE between the simulation temperature values and experimental temperature values as the fitness function, the CHTC optimization program was established. The population size was set as 20, maximum generation as 200, cross probability as 0.4, and mutation probability as 0.2; the optimal CHTCs were searched automatically without human intervention. The validation results showed that the simulation results agreed well with the experimental results with the maximum temperature simulation error decreased from 46.7% to 8.1% and axial (Z direction) thermal fluctuation simulation error decreased from 10.34% to −1.2%. The results demonstrated the efficiency of the proposed method for searching the optimal CHTCs especially those with many CHTCs. 4. From the aspect of optimizing CHTCs in FETA of the spindle system, the method proposed in this paper improved some accuracy, but there still exists a certain simulation error, which may be the uncertainty in motor heat and TCRs. So, next we may focus on optimizing the whole boundary conditions, i.e., heat generation, TCRs, and CHTCs, during FETA in the future.

3.

4.

5.

6.

7.

8. 9.

10.

11.

12.

13.

14. 15.

16.

17.

18.

Ramesh R, Mannan MA, Poo AN (2000) Error compensation in machine tools—a review: part II: thermal errors. Int J Mach Tools Manuf 40(9):1257–1284 Mayr J, Jedrzejewski J, Uhlmann E, Alkan Donmez M, Knapp W, Hartig F, Wendt K, Moriwaki T, Shore P, Schmitt R (2012) Thermal issues in machine tools. CIRP Ann Manuf Technol 61(2):771–791 Liang RJ, Ye WH, Zhang HH, Yang QF (2012) The thermal error optimization models for CNC machine tools. Int J Adv Manuf Technol 63(9–12):1167–1176 Wang BM, Mei XS, ZX W, Zhu F (2015) Dynamic modeling for thermal error in motorized spindles. Int J Adv Manuf Technol 78(5– 8):1141–1146 Creighton E, Honegger A, Tulsian A, Mukhopadhyaye D (2010) Analysis of thermal errors in a high-speed micro-milling spindle. Int J Mach Tools Manuf 50(4):386–393 ZZ X, Liu XJ, Kim HK, Shin JH, Lyu SK (2011) Thermal error forecast and performance evaluation for an air-cooling ball screw system. Int J Mach Tools Manuf 51(7):605–611 Liu ZF, Pan MH, Zhang AP, Zhao YS, Yang Y, Ma CY (2015) Thermal characteristic analysis of high-speed motorized spindle system based on thermal contact resistance and thermal-conduction resistance. Int J Adv Manuf Technol 76(9–12):1913–1926 Harris TA, Kotzalas MN (2006) Essential concepts of bearing technology. CRC press Sun LJ, Ren MJ, Hong HB, Yin YH (2016) Thermal error reduction based on thermodynamics structure optimization method for an ultra-precision machine tool. Int J Adv Manuf Technol:1–11 Burghold EM, Frekers Y, Kneer R (2015) Determination of timedependent thermal contact conductance through IR-thermography. Int J Therm Sci 98:148–155 Zhao HL, Huang YM, Xu JL, Jiang LY, Zhang WH, Sheng BH (1999) Experiment research on thermal contact resistance of normal used joints. J Xi’an Univ Technol 15(3):26–29 (in Chinese) Xu RP, Feng HD, Zhao LP, Xu L (2006) Experimental investigation of thermal contact conductance at low temperature based on fractal description. Int Commun Heat and Mass 33(7):811–818 Ma C, Mei XS, Yang J, Zhao L, Shi H (2015) Thermal characteristics analysis and experimental study on the high-speed spindle system. Int J Adv Manuf Technol 79(1–4):469–489 Holman JP (2010) Heat transfer. McGraw-Hill Min X, Jiang S (2011) A thermal model of a ball screw feed drive system for a machine tool. Proc Inst Mech Eng Part C J Mech Eng Sci 225(1):186–193 Yang H, Yin GF, Fang H, Liu LX, Zhang DJ, Xu DW (2011) Research on calculation methods of convective heat transfer coefficients for machine tools’ thermal analysis. J Sichuan Univ Eng Sci Edit 43(4):241–248 (in Chinese) Chow JH, Zhong ZW, Lin W, Khoo LP (2012) A study of thermal deformation in the carriage of a permanent magnet direct drive linear motor stage. Appl Therm Eng 48:89–96 Li DX, Feng PF, Zhang JF, Wu ZJ, Yu DW (2016) Method for modifying convective heat transfer coefficients used in the thermal simulation of a feed drive system based on the response surface methodology. Numer Heat Tr A-Appl 69(1):51–66

Int J Adv Manuf Technol 19.

20. 21.

22. 23.

ZC D, Yao SY, Yang JG (2015) Thermal behavior analysis and thermal error compensation for motorized spindle of machine tools. Int J Precis Eng Manuf 16(7):1571–1581 Xu M, Jiang SY, Cai Y (2007) An improved thermal model for machine tool bearings. Int J Mach Tools Manuf 47(1):53–62 Yang ZW, Yin GF, Shang X, Jiang H, Zhong KY (2011) Coupling analysis model of thermal and dynamic characteristics for highspeed motorized spindle. J Jilin Univ Eng Technol Edit 41(1): 100–105 (in Chinese) Bossmanns B, JF T (1999) A thermal model for high speed motorized spindles. Int J Mach Tools Manuf 39(9):1345–1366 Sardiñas RQ, Santana MR, Brindis EA (2006) Genetic algorithmbased multi-objective optimization of cutting parameters in turning processes. Eng Appl Artif Intell 19(2):127–133

24.

Gosselin L, Tye-Gingras M, Mathieu-Potvin F (2009) Review of utilization of genetic algorithms in heat transfer problems. Int J Heat Mass Tran 52(9):2169–2188 25. Wang GL, Wang YQ, Zhao J, Chen GL (2012) Process optimization of the serial-parallel hybrid polishing machine tool based on artificial neural network and genetic algorithm. J Intell Manuf 23(3):365–374 26. Wang JS, Zhu CA, Feng MC, Ren WQ (2013) Thermal error modeling and compensation of long-travel nanopositioning stage. Int J Adv Manuf Technol 65(1–4):443–450 27. Huang YQ, Zhang J, Li X, Tian LJ (2014) Thermal error modeling by integrating GA and BP algorithms for the high-speed spindle. Int J Adv Manuf Technol 71(9–12):1669–1675