Next Article in Journal
A Method for Dynamically Selecting the Best Frequency Hopping Technique in Industrial Wireless Sensor Network Applications
Next Article in Special Issue
Glucose Sensor Using U-Shaped Optical Fiber Probe with Gold Nanoparticles and Glucose Oxidase
Previous Article in Journal
Microwave Sensors for Breast Cancer Detection
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Parameter Estimation of the Thermal Network Model of a Machine Tool Spindle by Self-made Bluetooth Temperature Sensor Module

1
Institute of Applied Mechanics, National Taiwan University, No. 1, Sec. 4, Roosevelt Rd., Taipei 10617, Taiwan
2
Department of Mechanical and Electromechanical Engineering, National ILan University, No.1, Sec. 1, Shennong Rd., Yilan City, Yilan County 260, Taiwan
*
Author to whom correspondence should be addressed.
Sensors 2018, 18(2), 656; https://doi.org/10.3390/s18020656
Submission received: 24 January 2018 / Revised: 15 February 2018 / Accepted: 21 February 2018 / Published: 23 February 2018
(This article belongs to the Special Issue Selected Sensor Related Papers from ICI2017)

Abstract

:
Thermal characteristic analysis is essential for machine tool spindles because sudden failures may occur due to unexpected thermal issue. This article presents a lumped-parameter Thermal Network Model (TNM) and its parameter estimation scheme, including hardware and software, in order to characterize both the steady-state and transient thermal behavior of machine tool spindles. For the hardware, the authors develop a Bluetooth Temperature Sensor Module (BTSM) which accompanying with three types of temperature-sensing probes (magnetic, screw, and probe). Its specification, through experimental test, achieves to the precision ±(0.1 + 0.0029|t|) °C, resolution 0.00489 °C, power consumption 7 mW, and size Ø40 mm × 27 mm. For the software, the heat transfer characteristics of the machine tool spindle correlative to rotating speed are derived based on the theory of heat transfer and empirical formula. The predictive TNM of spindles was developed by grey-box estimation and experimental results. Even under such complicated operating conditions as various speeds and different initial conditions, the experiments validate that the present modeling methodology provides a robust and reliable tool for the temperature prediction with normalized mean square error of 99.5% agreement, and the present approach is transferable to the other spindles with a similar structure. For realizing the edge computing in smart manufacturing, a reduced-order TNM is constructed by Model Order Reduction (MOR) technique and implemented into the real-time embedded system.

1. Introduction

The thermal characteristics of machine tool spindles play an important role in the development of high-precision and high-speed machining. Many literatures about the machining error analysis had indicated that a large amount of the total machining errors attributed to the thermal issue [1,2]. In recent years, diagnostics and prognostics of sudden failure of machine tool spindle have drawn considerable attention. Hence, the model-based method is the key technology for predicting the thermal behavior of spindle. Bossmanns and Tu et al. developed a comprehensive thermo-mechanical model to characterize the heat transfer and quantitative power flow of a high-speed motorized spindle based on Finite Difference Method (FDM) [3,4,5]. On the other hand, Finite Element Method (FEM) has the advantage of fully simulating the thermo-dynamic-mechanical behavior of the spindle. These FEM models not only simulate the temperature distribution and the thermal deformation of machine tool, but also take the bearing stiffness and thermal contact resistance into consideration [6,7]. Experimental and simulation evidence indicated that the temperature distribution of spindle is axisymmetric [8]. Meanwhile, the thermal resistance network model was proposed considering the thermal contact resistance based on fractal and Hertz contact theory [9,10]. Several relevant issues were widely investigated in thermal characteristics of machine tool spindle, such as thermal contact resistance, lubricant viscosity, cooling condition, thermally induced preload, and transient thermal behavior [6,11,12,13]. Heat flux is one of the most important parameters of the thermal behavior of spindle. The heat flux of a thermal system could be inversely predicted by recursive least squares algorithm from the thermal resistance network model [14] or finite element model [15]. System identification techniques have drawn much attention to many thermal predictive and control applications, such as the core temperature real-time estimation for lithium ion batteries [16], the temperature prediction of the squirrel-cage rotor of induction motor [17], the thermal system analysis of buildings [18], and the transient thermal impedance of semiconductor [19]. Regarding the failure prevention of permanent magnet synchronous motors, the grey-box system identification of the model-based lumped-parameter thermal networks (LPTNs) has been developed to estimate the permanent magnet and winding temperature [20,21,22].
This article presents a simplified and reduced-order thermal model for a machine tool spindle, which can predict the thermal response of certain crucial points within the spindle. Though FEM and FDM are great approaches to analyze the thermal characteristic of the spindle; however, they require a large amount of numerical computation. These methods are unsuitable for in-situ application due to their high demand of numerical computation. Instead of comprehensively analyzing the heat transfer equations, this article establishes a lumped-element Thermal Network Model (TNM) that effectively reduces the computation load while retaining the predictive performance. The TNM consists of the parameters of thermal resistances, thermal capacitances, and heat sources. The assumptions and simplifications of those parameters are made based on the theory of heat transfer and empirical formula. Note that the global and local solutions of the parameters are solved from the measured temperature data by means of nonlinear least square method and grey-box system identification method. The TNM is verified by self-validation as well as external-validation. Moreover, the Short-Circuit Time Constant (SCTC) method is applied to approximate the low-frequency band limitation. Also, the Model Order Reduction (MOR) techniques allow for reducing the order of the estimated model while retaining the performance. Furthermore, Bluetooth Temperature Sensor Module (BTSM) is developed that provides high accuracy, high stability, and five-channel temperature measurement. Three types of temperature sensing probes (magnetic, screw, and probe) based on Resistance Temperature Detector (RTD) are constructed in this work, those are suitable for being placed at some characteristic locations such as into the outer-ring of bearings, approaching to the spindle, and on the surfaces of the housing.

2. Theoretical Background

2.1. Thermal Characteristics of Spindle

Many researchers have comprehensively investigated the thermal characteristics of machine tool spindle, including heat generation, forced heat convection over rotating shaft and annulus, heat conduction, heat radiation, and heat capacity. The primary concern is the speed-related influence, hence, the degree of each leading term related to the rotating speed should be determined especially. The frictional torque of the angular contact ball bearings causes the major heat generation of the externally driven spindle. Harris proposed that the heat sources were mainly caused by three factors: applied load (Ql), viscous shear (Qv), and spinning motion (Qs) [23].
Q f = Q l + Q v + Q s = 1.047 × 10 4 × [ n ( M l + M v ) + j = 1 Z n s i M s i , j + n s o M s o , j ]
where Ml is the applied-load-dependent term, Mv is the viscous-shear-dependent term, and Msi and Mso are the spinning-motion-dependent terms, n denotes the rotational speed of bearing, nsi and nso denotes the spin speeds of the inner and outer rolling element, respectively, and Z is the number of rolling elements.
Stein and Tu modified Palmgren’s equation for specifying the thermally induced preload of bearing [24]. The bearing friction torque due to applied load (Ml) is obtained by Equation (2). In their research, the predicted thermally induced preload would reach steady over time and its magnitude grew with the rotational speed. The nonlinear relationship between thermally induced preload and rotational speed is like the power function.
M l = μ l f l ( F n + F t ) ( r i + r b )
where Ft is the thermally induced preload, Fn is the contact force from external load, ri and rb is the radius of inner ring and rolling element, μl is the friction coefficient due to load. The empirical Equation (3) explains the frictional torque that is caused by lubricant shear viscosity (Mv) within the bearing.
M v = 10 7 f 0 ( ν L n ) 2 / 3 d m 3 ,   for   ν L n 2000 M v = 160 × 10 7 f 0 d m 3 ,   for   ν L n < 2000
where f0 is a factor depending on bearing type and lubrication, νL is lubricant kinematic viscosity. For angular contact ball bearing, the gyroscopic moment of the rolling elements necessarily leads to boring motion that would cause the friction torque on the contact surface, and thus also introduce heat generation. The friction torque related to spinning motion (Ms) is obtained by Equation (4) [25].
M s = 3 μ s F r o l l a e 8
where Froll is the contact load, a represent the major axes of the elliptical contact area, e is the elliptical integral of the second kind, and μs is the friction coefficient.
The primary concern is to evaluate the heat source within the spindle thermal system associated with the spindle rotational speed (ω). Many researchers had indicated that the viscous friction has the influential contribution in total frictional torque within the bearing [12]. As for simplification, due to the power function trend of thermally induced preload influence, the polynomial-like function might have the ability to predict the heat generation of bearing. As mentioned above, substituting Equations (2)–(4) into Equation (1), the authors conclude that the heat source can be expressed as a function of the viscous-related term (ω5/3) and the remaining term (ω), as described in Equation (5).
Q f ~ f ( ω , F t )      ( 1 )      f ( ω 5 / 3 ,   ω 1 ) ( 1 ) F t ~ f ( ω c 0 )
An analytical solution describes the forced convective heat transfer over rotating shaft [26]. The average Nusselt number ( Nu ¯ ), Reynolds number (Re), and Prandtl number (Pr) are determined by the following equations.
Nu ¯ = 0.6366 ( Re Pr ) 1 / 2
Re = ω d 2 2 ν a i r , Pr = ν a i r α a i r
Assume constant material properties, the rotational speed dominates the forced heat convection. This implies that the forced convection coefficient is a speed-related term being proportional to ω0.5. In this case, we take the rotational speed as the leading term and express the thermal resistance of forced convection around moving surface as a function of ω−0.5.
h f o r c e d = Nu ¯ × k a i r d Re 0.5 ω 0.5
R f o r c e d = 1 h f o r c e d A ~ f ( ω 0.5 )
The forced heat convection over the surfaces of rotating annulus can be determined by the following Equations (10) [27]. Taylor number (Ta) and the geometry factor (Fg, P) are obtained. Note that the geometry factors (Fg and P) depend on the inner and outer radius of annulus only; therefore, it can be viewed as a constant value.
Nu ¯ = 0.409 ( Ta 2 F g ) 0.241 ,   for   10 4 < Ta 2 F g < 10 7
where
Ta = ω [ ( r o + r i ) 2 ] 0.5 ( r o r i ) 1.5 ν ; F g = ( π 4 1697 ) ( 1 ( r o r i ) ( r o + r i ) ) 2 P 1 ;
P = 0.0571 { 1 0.652 ( r o r i r i ) } + 0.00056 { 1 0.652 ( r o r i r i ) } 1
and ri and ro are the inner and outer radius. Namely, the coefficient of forced heat convection can be derived in Equation (11), which is proportional to ω0.482. Thus, it is appropriate to formulate the thermal resistance of forced convection near annulus as a function of ω−0.482 in Equation (12).
h a n n u l u s Nu ¯ Ta 0.482 ω 0.482
R a n n u l u s = 1 h a n n u l u s A a n n u l u s ω 0.482
According to [28], free convection around a horizontal cylinder and stationary ambient can be estimated by the following Equation (13).
Nu ¯ = { 0.6 + 0.387 [ Gr Pr [ 1 + ( 0.559 / Pr ) 9 / 16 ] 16 / 9 ] 1 / 6 } 2   f o r   10 5 < Gr Pr < 10 12
where
Pr = ν a i r α a i r ; Gr = g β ( T s T a ) D 3 ν a i r 2 = ( 1 ) g ( T s T a ) D 3 T f ν a i r 2 ( 1 ) β ~ 1 T f ( i d e a g a s a p p r o x i m a t i o n ) ; T f = T s + T a 2
Except for the material properties and geometry factor, the free convection coefficient is dependent on the spindle surface temperature (Ts) and ambient temperature (Ta) particularly. It would be difficult to simulate the free convective phenomenon according to the temperature difference. Alternatively, the rotational speed has strong correlation to the temperature difference. We assume that the free convective coefficient can be expressed as the function in Equation (14). Hence, the free convective thermal resistance can be derived as two terms, a constant term and a speed-related term. Especially, all of the assumptions must be verified by experimental results, and it will be further discussed in Section 5.
h f r e e = h f r e e , c + h f r e e , s ω n R f r e e = 1 h f r e e , c A + h f r e e , s ω n A = R f r e e , c R f r e e , s
In general, the radiative thermal resistance can be expressed as Equation (15). It indicates that the resistance of thermal radiation depends on surface temperature (Ts) and ambient temperature (Ta). Due to the extremely tiny value of Stefan-Boltzmann constant (σ) and the emissivity (ε), the influence of surface temperature variation would be limited. Consequently, the radiative thermal resistance can be considered as a constant value, which would be further verified by experiment.
R r a d = 1 ε σ A ( T s 2 + T a 2 ) ( T s + T a )
The thermal capacitance is directly determined by the element properties (volume, heat capacity, density) through Equation (16). It represents the store of energy that is required to increase the temperature. In other words, the thermal capacitance value is determined once the spindle structure is designed in practice. Thereby, it is assumed as a constant, when considering that the material properties remain relatively the same under different operating conditions. The estimated value of each thermal capacitance has the decisive influence on the thermal time constant.
C = ρ V c p

2.2. System Identification Technique

The TNM of machine tool spindle is identified by system identification methodology based on physical concept, namely the so-called grey-box estimation. It is utilized to predict the temperature distribution of spindle. Figure 1 explains the system identification procedure for the TNM of machine tool spindle comprehensively. First, through the temperature data measured at steady state, solve for the thermal resistances of TNM at operating mode by means of the “lsqnonlin” function and the “trust-region-reflective” algorithm; the aforesaid function and algorithm are the nonlinear least square regression provided by MATLAB [29]. Second, solve for the thermal capacitance of TNM by means of the “greyest” function provided by MATLAB; that function is exactly the grey-box model identification [29]. Finally, the thermal resistances at natural cooling mode are obtained by the grey-box model identification as well. The least square system identification problem usually exists non-unique solutions satisfying the cost function (minimization). It is expected that some irrational solutions can be eliminated based on physical sense, for example, the convective thermal resistance is generally greater than the conductive one. Furthermore, it must be emphasized that while inversely estimate the heat generation, thermal resistance, and thermal capacitance at the same time, some analogous solutions might occur because different solutions with the same ratio among all of the parameters might introduce the same result. Fortunately, the total thermal capacitance ( j = 1 4 C j ) is a constant value; thus, the appropriate solution can be determined accordingly. Finally, it is necessary to verify the performance of the estimated TNM; hence, the self-validation and external validation are further investigated in Section 5.

3. Parameterization Methodology for the Thermal Network Model of Spindle

3.1. Parameterization Strategy

The assumptions for modeling the TNM are listed as follows:
(1)
lumped element model assumption is valid, Bi << 0.1;
(2)
temperature distribution of spindle is axisymmetric;
(3)
heat generation and forced convective resistance are assumed according to the theoretical and empirical Equations (5), (9), and (12);
(4)
heat transfer through radiation and conduction are considered as a constant thermal resistance, including thermal contact resistance; and,
(5)
free convection coefficient is assumed as a function of h f r e e , c + h f r e e , s ω 0.5 .
According to previous research, the temperature distribution of spindle is almost axisymmetric. For this reason, the two-dimensional heat equation with internal heat generation can be implemented to describe the three-dimensional thermal behavior. Furthermore, based on the assumption of lumped element model, the TNM can characterize the heat transfer inside the spindle.

3.2. Estimated Thermal Network Model in State-Space

While the spindle operates at a constant rotational speed, the frictional heat generation are assumed to be a constant source, denoted as Qf1, Qf2. As a result, the heat balance equation with internal heat generation is expressed as
[ C ] { T } t = [ K ] { T } + [ K a ] { T a } + { Q H }
The authors define the temperature difference (θ) with respect to the ambient temperature.
{ θ ( t ) } { T ( t ) } T a t m
Substituting it into Equation (17) gives the heat balance equations
[ C ] { θ } t = [ K ] { θ } + { Q H }
Consider the physical insight and geometry of the spindle; we construct the lumped-parameter TNM of the spindle with four nodes, illustrated in Figure 2, which are rear bearing A (T1), midpoint of housing (T2), front bearing D (T3), midpoint of shaft (T4), respectively. The specific locations will be illustrated latterly in Section 4.2. TNM contains four types of components, thermal resistance, speed-related thermal resistance, thermal capacitance, and heat source. The authors classify the thermal behavior of spindle into two modes, operating mode and natural cooling mode. At operating mode, the spindle is operated at a constant rotation speed and the initial temperature is equal to the ambient temperature, illustrated in Figure 2a. At natural cooling mode, the spindle is stopped and releases the heat that is stored during operating mode, as illustrated in Figure 2b. The state-space representation of the TNM at operating mode is expressed as
{ x ˙ ( t ) = A o p x ( t ) + B o p u ( t ) y ( t ) = C o p x ( t ) + D o p u ( t )
According to the assumptions, the convective resistance coefficient (ravj, r24) and the heat generation coefficient (qf1qf4) are unknown parameters needed to be solved.
{ R a v j = r a v j ( ω / ω 0 ) 0.5 , R 24 = r 24 ( ω / ω 0 ) 0.482 Q f 1 = q f 1 ( ω / ω 0 ) 5 / 3 + q f 3 ( ω / ω 0 ) 1 Q f 2 = q f 2 ( ω / ω 0 ) 5 / 3 + q f 4 ( ω / ω 0 ) 1
where ω0 is 5000 rpm, thus the unit of the coefficients can be meaningful. At steady state, the thermal capacitances can be ignored due to the temperature derivative terms decay at steady state. Therefore, rewriting the heat balance equation as:
A s t d θ + B s t d θ ω 0.482 + C s t d θ ω 0.5 + D s t d ω 5 / 3 + E s t d ω 1 = 0
Equation (22) is a Single-Input-Multiple-Output (SIMO) system whose input is the spindle rotation speed (ω) while the spindle temperatures (θ) are the outputs. As the speed-related parameters varying with different spindle rotation speeds, the thermal response of TNM would change as well; therefore, these nonlinear behaviors make the TNM more complicated. Likewise, the state-space representation of the TNM at natural cooling mode is expressed as:
{ x ˙ ( t ) = A n c x ( t ) + B n c u ( t ) , x i ( 0 ) = θ i ( 0 ) y ( t ) = C n c x ( t ) + D n c u ( t )
The specific system parameter matrices are shown in Appendix.

4. Experiment Setup

4.1. Self-Made Bluetooth Temperature Sensor Module

The performance (accuracy, resolution, stability…) of the sensor is one of the key points directly affecting the quality of the system identification results. Due to the high stability and accuracy of RTD, the authors develop a BTSM that provides multi-channel temperature measurement with high precision. Figure 3a shows the circuit design of BTSM, it is implemented with a Bluno Beetle [30] as the microcontroller, a 24-bit AD7794 [31] as the analog-to-digital converter, a reference resistor with low temperature coefficient of resistance, 4-wire Kelvin measurement, and the Printed Circuit Board (PCB) Layout (Figure 3b). After several versions of evolution, the BTSM (Figure 4a) has achieved with high accuracy of ±(0.1 + 0.0029|ϑ|) °C, the resolution of 0.00489 °C, and the module size is only Ø40 mm × 27 mm. The specification is listed in Table 1. As the front-end temperature acquisition element, three types of temperature probe (screw type, magnetic mount type, and probe type) are fabricated by platinum RTD (PT1000) and sealed with thermal conductive Gel (Figure 4a). Magnet type probe is made for surface temperature measurement, and the other types are made for the shaft, housing, bearing temperature measurement. The fabrication process of the temperature probe is shown in supplementary materials Figure S1. It should be mentioned here that, to improve the precision of the BTSM, the Kalman filter is used to filter the noise from the measured temperature signal.

4.2. Experiment Setup

For validating the commercial application of the present TNM and BTSM, the authors purchased a customized, but without loss of generality, spindle warm-up system. Illustrated in Figure 4b, it contains a spindle driven by a 7.5 kW motor through belt, a controller accompanying with a PC being used to conduct the programmable control of motor, a commercial thermocouple device being used to verify the feasibility of the self-made BTSM, etc. The embedded locations of the temperature probes are designed based on the thermal characteristic of machine tool spindle. Specifically illustrated in Figure 5, the temperature sensors are strategically located at 12 nodes near the front and rear bearings, spindle housing, shaft, spindle surface, and ambient temperature, respectively.

5. Results

The transient temperature distribution of spindle is measured by BTSM at various rotation speeds. Figure 6a shows a typical temperature rise curves at the rotation speed of 6021 rpm. Note that, at operating mode is similar to the charging mode of R-C circuit, more specifically, the zero-state response. Likewise, the temperature fall behavior at the natural cooling mode is similar to the discharging mode of R-C circuit, more specifically, the zero-input response. Turn off the spindle right after reaching its steady state, meanwhile insert the temperature probe into the shaft. Consequently, the complete steady-state temperature distribution of spindle can be measured. Zoom in the region of the rectangular frame line of Figure 6a, namely the transition from operating mode to natural cooling mode, Figure 6b reveals that the uneven temperature distribution accumulated at operating mode will gradually reach the thermal equilibrium in nearly 1000 s. After that, the entire spindle behaves relatively as a homogeneous temperature element and cools through free heat convection and radiation heat transfer. It implies that the thermal convective resistances are greater than the thermal conductive resistances. As a result, it further proves that the lumped element assumption is valid. Following the system identification procedure, as mentioned in Section 2.2, the estimated parameters of the TNM are listed in Table 2 and Table 3, and the value of the thermal parameters are calculated at the rotational speed of 5000 rpm.
To verify the assumptions of the free convective and radiative heat transfer, the measured temperature of rear bearing (T1) at steady state is used to simulate hfree and hrad. Due to the temperature of rear bearing is the highest one among all of the other positions, the simulation results can indicate the significant influence. In Figure 7a, the approximate computation of Rrad (the purple dashed line) from Equation (15) is almost a constant value. As for free convective heat transfer coefficient, the authors consider the spindle to be a homogeneous cylinder with the diameter (D) of 0.135 m and length (L) of 0.278 m, and its other material properties are based on 25 °C. Consequently, the approximately calculated Rfree (green dotted line) from Equation (13) can be curve-fitted with the formula 1 / A ( h f r e e , c + h f r e e , s ω 0.5 ) . For this reason, it is appropriate to adopt these assumptions. Figure 7b shows the predicted heat generations of rear and front bearings. One can estimate the dissipated heat through radiation (Qrad.approx.) and free convection (Qfree.approx.) by means of the preceding approximated hfree and hrad. Furthermore, Figure 7a shows the speed-dependent parameters with respect to the rotational speed. It indicates that the nonlinear behavior is more significant under lower rotational speed. For this reason, it is necessary to take the speed-dependent thermal behavior into consideration. In other words, due to the varying of the thermal parameters, the poles of the estimated TNM would slightly change with different operating speeds, which will be discussed latterly in Section 5.4.

5.1. Steady State Self-Validation

After applying the nonlinear least square estimation to find the minimum 2-norm solution of the overdetermined system in state space representation, the global solutions at steady state can be found, which represent the sets of the optimized thermal resistance value. First, the model is self-validated at steady state through a mutual comparison of the predicted and measured temperature at five different rotation speeds, say 4006, 5013, 6021, 7023, and 8028 rpm, as shown in Figure 8a. The predicted and measured temperatures agree very well with each other, and the maximum mean average percentage error (MAPE) is only 5.54%. Moreover, it should be emphasized that the heat generation and steady state temperature are closely related, which is confirmed from the similarity of Figure 7b and Figure 8b.

5.2. Transient State Self-Validation

The grey-box model identification technique is implemented to estimate the remaining parameters of the TNM at both operating and natural cooling modes. The transient self-validations are obtained by comparing the transient temperatures rise estimated by TNM with the measured ones. Five transient states are validated; namely, the spindle runs from stationary to five different steady speeds, say 4006, 5013, 6021, 7023, and 8028 rpm. Figure 9 and Figure 10 show promising validation with the minimum normalized mean square errors (NMSE) of 94.2% at operating mode and 95.8% at natural cooling mode. When the spindle operates at operating mode, the spindle shaft rotates at a constant speed; hence, the temperature probe is not able to contact the shaft. For this reason, the experimental temperature of spindle shaft (T4) cannot be obtained, and the T4-est represents the estimated shaft temperature, as shown in Figure 9d.

5.3. External Validation

To verify the performance and robustness of the TNM, the external validation is necessary. The experiment under stepwise rotational speed (3001, 5018, 7028 rpm) is utilized especially. The results (Figure 11) reveal that the predicted temperature is in satisfactory agreement with the experimental results. As a result, this TNM methodology is appropriate to implement to the machine tool spindle as a simplified thermal model.

5.4. Model Order Reduction

To implement the system parameter estimation into a micro control unit for real-time onboard application on edge computing, the order of TNM model should be further simplified to reduce the load of computation and data storage. MOR technique is an excellent approach for reducing the complexity of the estimated TNM. Balanced realization method is proposed to reduce the 4th-order TNM into 1st-order truncated model, which is a simplified realization with equal controllability and observability. Figure 12a shows the Hankel singular values of the estimated TNM; it points out that there exists a leading state with the larger value than other states. For this reason, we remove the less three states (2, 3, 4) and apply the first dominant state to reconstruct the 1st-order truncated TNM. The pole location of the estimated TNM and truncated model are shown in Figure 12c. In addition, the Bode diagram illustrates the relationship between the input heat generation (Qf1, Qf2) and the output temperature T1 to T4 at rotational speed of 6000 rpm. It indicates that the truncated model has equivalent performance to the original estimated TNM under low frequency (<10−3) condition, as shown in Figure 13.

5.5. Short Circuit Time Constant

The Short Circuit Time Constant (SCTC) method provides an approximation of the poles of the estimated TNM. Each pole can be calculated according to each thermal capacitance multiplied by the nearby equivalent parallel thermal resistance while the other capacitances are short circuited, refer to Figure 12b. As the aspect of system degradation, we can simply recognize the time-varying trend by examining the variation of the time constant instantly by applying SCTC method.

6. Conclusions

In this article, we investigate the possibility of simplifying the spindle thermal behaviors by estimating the lumped-parameter TNM with system identification technique. We classify the spindle heat transfer condition into two modes, the operating mode and the natural cooling mode. The predicted thermal response of TNM agrees well with the experimental temperature results, verified by self-validation and external validation. In the aspect of the real-time onboard application, the reduced-order model by balanced realization has successfully lowered the computational complexity. In practice, the spindle warm-up procedure is essential to all of the machine tool spindles, especially for the new spindle. By operating at step-wise speeds, the lubricant within the spindle bearings can be evenly spread. The modeling methodology in this research can be appropriately integrated with the warm-up procedure and establish the estimated TNM for further operation. Figure 14 summarizes that the TNM methodology has the ability to predict the spindle thermal behavior; and the BTSM provides real-time monitoring the spindle temperature distribution. In that case, the vision of edge computing within a minimize temperature sensor module can be realized by integrated with the estimated TNM performed as the thermo-feature identification system (TID-sys), which can be customized for each machine tool spindle. Figure 15 demonstrates the real-time user interface of TID-sys on smartphone and PC. On the right-hand side, the TID-sys diagnosis the state (normal or alarm) of the machine tool spindle displayed through the terminal App [32] on smartphone. On the left-hand side, the measured temperature (black dot) of rear bearing (T1), midpoint of housing (T2), and front bearing (T3) attempt to grow along with the curve predicted by the estimated TNM. The TID-sys ensure that the machine tool spindle runs under normal condition. Once the significant temperature variation or abnormal temperature rise is detected, the operator can slow down the spindle or even stop for early protection. The full demonstration video is shown in supplementary materials Video S1.

Supplementary Materials

The following are available online at https://www.mdpi.com/1424-8220/18/2/656/s1, Figure S1: Fabrication process of the temperature probe, Figure S2: Experimental setup; and Video S1: Demonstration of the TID-sys.

Acknowledgments

This work was supported by the Ministry of Science and Technology, Taiwan through the grant number: MOST 105-2221-E-197-012 and MOST 106-2221-E-197-014.

Author Contributions

This paper was performed in collaboration among the authors. Yuan-Chieh Lo was responsible for the parameter estimation of the thermal network model, experiment works, and writing the article. Yuh-Chung Hu was the research project instructor and responsible for advising the research works and article writing. Pei-Zen Chang is Lo’s graduate advisor.

Conflicts of Interest

The authors declare no conflict of interest.

Appendix A. System Parameter Matrices

The input, output and state vector of TNM and linear in state space representation:
x = [ θ 1 θ 2 θ 3 θ 4 ] T u = [ Q f 1 0 Q f 2 0 ] T y = [ θ 1 θ 2 θ 3 θ 4 ] T
The heat balance equations of TNM:
{ θ 1 R a 1 + θ 1 R a v 1 + θ 2 θ 1 R 12 + θ 4 θ 1 R 14 + Q f 1 = C 1 θ ˙ 1 θ 2 R a 2 + θ 1 θ 2 R 12 + θ 3 θ 2 R 23 + θ 4 θ 2 R 24 = C 2 θ ˙ 2 θ 3 R a 3 + θ 3 R a v 2 + θ 2 θ 3 R 23 + θ 4 θ 3 R 34 + Q f 2 = C 3 θ ˙ 3 θ 4 R a v 3 + θ 1 θ 4 R 14 + θ 2 θ 4 R 24 + θ 3 θ 4 R 34 = C 4 θ ˙ 4
The TNM system matrices of state-space representation at operating mode:
A o p = A o p _ c o n d + A o p _ c o n v A o p _ c o n d = [ ( R 12 1 R 14 1 ) C 1 1 R 12 1 C 1 1 0 R 14 1 C 1 1 R 12 1 C 2 1 ( R 12 1 R 23 1 R 24 1 ) C 2 1 R 23 1 C 2 1 R 24 1 C 2 1 0 R 23 1 C 3 1 ( R 23 1 R 34 1 ) C 3 1 R 34 1 C 3 1 R 14 1 C 4 1 R 24 1 C 4 1 R 34 1 C 4 1 ( R 14 1 R 24 1 R 34 1 ) C 4 1 ] A o p _ c o n v = [ ( R a 1 1 R a v 1 1 ) C 1 1 0 0 0 0 R a 2 1 C 2 1 0 0 0 0 ( R a 3 1 R a v 2 1 ) C 3 1 0 0 0 0 R a v 3 1 C 4 1 ] B o p = d i a g ( C 1 1 , C 2 1 , C 3 1 , C 4 1 ) ; C o p = I 4 × 4 ;   D o p = O 4 × 4 ;
The TNM system matrices of state-space representation at steady state:
A s t d = [ R a 1 1 R 12 1 R 14 1 R 12 1 0 R 14 1 R 12 1 R a 2 1 R 12 1 R 23 1 R 23 1 0 0 R 23 1 R a 3 1 R 23 1 R 34 1 R 34 1 R 14 1 0 R 34 1 R 14 1 R 34 1 ] B s t d = [ 0 0 0 0 0 r 24 1 0 r 24 1 0 0 0 0 0 r 24 1 0 r 24 1 ] ; C s t d = [ r a v 1 1 0 0 0 0 0 0 0 0 0 r a v 2 1 0 0 0 0 r a v 3 1 ] ; D s t d = [ q f 1 q f 1 0 0 q f 2 q f 2 0 0 ] ; E s t d = [ q f 3 q f 3 0 0 q f 4 q f 4 0 0 ] ;
The TNM system matrices of state-space representation at natural cooling mode:
A n c = A n c _ c o n d + A n c _ c o n v A n c _ c o n d = [ ( R 12 1 R 14 1 ) C 1 1 R 12 1 C 1 1 0 R 14 1 C 1 1 R 12 1 C 2 1 ( R 12 1 R 23 1 R 24 1 ) C 2 1 R 23 1 C 2 1 R 24 1 C 2 1 0 R 23 1 C 3 1 ( R 23 1 R 34 1 ) C 3 1 R 34 1 C 3 1 R 14 1 C 4 1 R 24 1 C 4 1 R 34 1 C 4 1 ( R 14 1 R 24 1 R 34 1 ) C 4 1 ] A n c _ c o n v = [ R a 1 1 C 1 1 0 0 0 0 R a 2 1 C 2 1 0 0 0 0 R a 3 1 C 3 1 0 0 0 0 R a 4 1 C 4 1 ] B n c = O 4 × 4 ;   C n c = I 4 × 4 ;   D n c = O 4 × 4 ;
Mean average percentage error:
MAPE = 1 n i = 1 n | x f o r e c a s t , i x a c t u a l , i x a c t u a l , i | × 100 %
Normalized mean square error:
NMSE = 1 x r e f x 2 x r e f x ¯ r e f 2

References

  1. Mayr, J.; Jedrzejewski, J.; Uhlmann, E.; Donmez, A.; Knapp, W.; Härtig, F.; Wendt, K.; Moriwaki, T.; Shore, P.; Schmitt, R.; et al. Thermal issues in machine tools. CIRP Ann. Manuf. Technol. 2012, 61, 771–791. [Google Scholar] [CrossRef]
  2. Ramesh; Mannan, M.; Poo, A. Error compensation in machine tools—A review. Part II: Thermal errors. Int. J. Mach. Tools Manuf. 2000, 40, 1257–1284. [Google Scholar] [CrossRef]
  3. Bossmanns, B.; Tu, J. A thermal model for high speed motorized spindles. Int. J. Mach. Tools Manuf. 1999, 39, 1345–1366. [Google Scholar] [CrossRef]
  4. Bossmanns, B.; Tu, J. A Power Flow Model for High Speed Motorized Spindles—Heat Generation Characterization. J. Manuf. Sci. Eng. 2001, 123, 494. [Google Scholar]
  5. Lin, C.-W.; Tu, J.; Kamman, J. An integrated thermo-mechanical-dynamic model to characterize motorized machine tool spindles during very high speed rotation. Int. J. Mach. Tools Manuf. 2003, 43, 1035–1050. [Google Scholar]
  6. Ma, C.; Yang, J.; Zhao, L.; Mei, X.; Shi, H. Simulation and experimental study on the thermally induced deformations of high-speed spindle system. Appl. Therm. Eng. 2015, 86, 251–268. [Google Scholar] [CrossRef]
  7. Zivkovic, A.; Zeljkovic, M.; Tabakovic, S.; Milojevic, Z. Mathematical modeling and experimental testing of high-speed spindle behavior. Int. J. Adv. Manuf. Technol. 2015, 77, 1071–1086. [Google Scholar] [CrossRef]
  8. Huang, Y.-H.; Huang, C.-W.; Chou, Y.-D.; Ho, C.-C.; Lee, M.-T. An Experimental and Numerical Study of the Thermal Issues of a High-speed Built-in Motor Spindle. Smart Sci. 2016, 4, 1–7. [Google Scholar] [CrossRef]
  9. Liu, Z.; Pan, M.; Zhang, A.; Zhao, Y.; Yang, Y.; Ma, C. Thermal characteristic analysis of high-speed motorized spindle system based on thermal contact resistance and thermal-conduction resistance. Int. J. Adv. Manuf. Technol. 2015, 76, 1913–1926. [Google Scholar]
  10. Brecher, C.; Shneor, Y.; Neus, S.; Bakarinow, K.; Fey, M. Thermal Behavior of Externally Driven Spindle: Experimental Study and Modelling. Engineering 2015, 7, 73–92. [Google Scholar]
  11. Min, X.; Shuyun, J.; Ying, C. An improved thermal model for machine tool bearings. Int. J. Mach. Tools Manuf. 2007, 47, 53–62. [Google Scholar] [CrossRef]
  12. Takabi, J.; Khonsari, M.M. Experimental testing and thermal analysis of ball bearings. Tribol. Int. 2013, 60, 93–103. [Google Scholar] [CrossRef]
  13. Wu, L.; Tan, Q. Thermal Characteristic Analysis and Experimental Study of a Spindle-Bearing System. Entropy 2016, 18, 271. [Google Scholar] [CrossRef]
  14. Noh, J.-H.; Cha, K.-U.; Ahn, S.-T.; Yook, S.-J. Prediction of time-varying heat flux along a hollow cylindrical tube wall using recursive input estimation algorithm and thermal resistance network method. Int. J. Heat Mass Trans. 2016, 97, 232–241. [Google Scholar] [CrossRef]
  15. Than, V.-T.; Wang, C.-C.; Ngo, T.-T.; Huang, J. Estimating time-varying heat sources in a high speed spindle based on two measurement temperatures. Int. J. Therm. Sci. 2017, 111, 50–65. [Google Scholar] [CrossRef]
  16. Lin, X.; Perez, H.; Siegel, J.; Stefanopoulou, A.; Li, Y.; Anderson, R.; Ding, Y.; Castanier, M. Online Parameterization of Lumped Thermal Dynamics in Cylindrical Lithium Ion Batteries for Core Temperature Estimation and Health Monitoring. IEEE Trans. Control Syst. Technol. 2013, 21, 1745–1755. [Google Scholar]
  17. Kral, C.; Habetler, T.; Harley, R.; Pirker, F.; Pascoli, G.; Oberguggenberger, H.; Fenz, C. Rotor Temperature Estimation of Squirrel-Cage Induction Motors by Means of a Combined Scheme of Parameter Estimation and a Thermal Equivalent Model. IEEE Trans. Ind. Appl. 2004, 40, 1049–1057. [Google Scholar] [CrossRef]
  18. Kim, D.; Cai, J.; Ariyur, K.; Braun, J. System Identification for Building Thermal Systems under the Presence of Unmeasured Disturbances in Closed Loop Operation: Lumped Disturbance Modeling Approach. Build. Environ. 2016, 107, 169–180. [Google Scholar] [CrossRef]
  19. Skibinski, G.L.; Sethares, W.A. Thermal parameter estimation using recursive identification. IEEE Trans. Power Electron. 1991, 6, 228–239. [Google Scholar] [CrossRef]
  20. Huber, T.; Peters, W.; Böcker, J. Monitoring critical temperatures in permanent magnet synchronous motors using low-order thermal models. In Proceedings of the International Power Electronics Conference (IPEC-Hiroshima 2014-ECCE ASIA), Hiroshima, Japan, 18–21 May 2014. [Google Scholar]
  21. Wallscheid, O.; Böcker, J. Global Identification of a Low-Order Lumped-Parameter Thermal Network for Permanent Magnet Synchronous Motors. IEEE Trans. Energy Convers. 2016, 31, 354–365. [Google Scholar] [CrossRef]
  22. Gaona, D.; Wallscheid, O.; Böcker, J. Glocal Identification Methods for Low-Order Lumped-Parameter Thermal Networks Used in Permanent Magnet Synchronous Motors. In Proceedings of the IEEE International Conference on Power Electronics and Drive Systems, Honolulu, HI, USA, 12–15 December 2017. [Google Scholar]
  23. Harris, T.A.; Kotzalas, M.N. Essential Concepts of Bearing Technology; Taylor & Francis: New York, NY, USA, 2006. [Google Scholar]
  24. Stein, J.L.; Tu, J.F. A state-space model for monitoring thermally induced preload in anti-friction spindle bearings of high-speed machine tools. ASME J. Dyn. Syst. Meas. Control 1994, 116, 372–386. [Google Scholar] [CrossRef]
  25. Harris, T.A.; Kotzalas, M.N. Advanced Concepts of Bearing Technology: Rolling Bearing Analysis; Taylor & Francis: New York, NY, USA, 2006. [Google Scholar]
  26. Kendoush, A. An Approximate Solution of the Convective Heat Transfer from an Isothermal Rotating Cylinder. Int. J. Heat Fluid Fl. 1996, 17, 439–441. [Google Scholar] [CrossRef]
  27. Childs, P.; Long, C. A Review of Forced Convective Heat Transfer in Stationary and Rotating Annuli. Proc. Inst. Mech. Eng. Part. C J. Mech. Eng. Sci. 1996, 210, 123–134. [Google Scholar] [CrossRef]
  28. Kothandaraman, C.P. Fundamentals of Heat and Mass Transfer; John Wiley & Sons: New York, NY, USA, 2006. [Google Scholar]
  29. Ljung, L. System Identification Toolbox—User’s Guide; MathWorks: Natick, MA, USA, 1999. [Google Scholar]
  30. Bluno Beetle SKU:DFR0339. Available online: https://www.dfrobot.com/wiki/index.php/Bluno_Beetle_SKU:DFR0339 (accessed on 12 December 2017).
  31. AD7794, 6-Channel, Low Noise, Low Power, 24-Bit Sigma Delta ADC. Available online: http://www.analog.com/en/products/analog-to-digital-converters/ad7794.html#product-overview (accessed on 12 December 2017).
  32. Bluno Terminal. Available online: https://itunes.apple.com/tw/app/bluno-terminal/id794109935?mt=8 (accessed on 12 December 2017).
Figure 1. System identification procedure.
Figure 1. System identification procedure.
Sensors 18 00656 g001
Figure 2. (a) 4-node Thermal Network Model (TNM) at operating mode; (b) 4-node TNM at natural cooling mode.
Figure 2. (a) 4-node Thermal Network Model (TNM) at operating mode; (b) 4-node TNM at natural cooling mode.
Sensors 18 00656 g002
Figure 3. (a) Circuit diagram of Bluetooth Temperature Sensor Module (BTSM); (b) Printed Circuit Board (PCB) layout of BTSM.
Figure 3. (a) Circuit diagram of Bluetooth Temperature Sensor Module (BTSM); (b) Printed Circuit Board (PCB) layout of BTSM.
Sensors 18 00656 g003
Figure 4. (a) Evolution of BTSM and three types of temperature probes; (b) Experimental setup.
Figure 4. (a) Evolution of BTSM and three types of temperature probes; (b) Experimental setup.
Sensors 18 00656 g004
Figure 5. Sensor location and TNM representation in machine tool spindle.
Figure 5. Sensor location and TNM representation in machine tool spindle.
Sensors 18 00656 g005
Figure 6. (a) Experimental result of 12 nodes on 6021 rpm; (b) The transition region of reaching thermal equilibrium.
Figure 6. (a) Experimental result of 12 nodes on 6021 rpm; (b) The transition region of reaching thermal equilibrium.
Sensors 18 00656 g006
Figure 7. (a) Simulated heat generation varying with rotational speeds; (b) Estimated convective resistance varying with rotational speeds.
Figure 7. (a) Simulated heat generation varying with rotational speeds; (b) Estimated convective resistance varying with rotational speeds.
Sensors 18 00656 g007
Figure 8. (a) Self-validation between predicted and measured steady state temperature, T1 to T4 are the positions of the rear bearing A, midpoint of inner housing, front bearing D, and midpoint of shaft, respectively; (b) Steady-state temperature validation.
Figure 8. (a) Self-validation between predicted and measured steady state temperature, T1 to T4 are the positions of the rear bearing A, midpoint of inner housing, front bearing D, and midpoint of shaft, respectively; (b) Steady-state temperature validation.
Sensors 18 00656 g008
Figure 9. Model self-validation at operating mode with various rotation speeds (a) Rear bearing (T1) temperature validation (b) Midpoint of housing (T2) temperature validation (c) Front bearing (T3) temperature validation (d) Midpoint of shaft (T4) temperature validation.
Figure 9. Model self-validation at operating mode with various rotation speeds (a) Rear bearing (T1) temperature validation (b) Midpoint of housing (T2) temperature validation (c) Front bearing (T3) temperature validation (d) Midpoint of shaft (T4) temperature validation.
Sensors 18 00656 g009
Figure 10. Model self-validation at natural cooling mode with various rotation speeds (a) Rear bearing (T1) temperature validation (b) Midpoint of housing (T2) temperature validation (c) Front bearing (T3) temperature validation (d) Midpoint of shaft (T4) temperature validation.
Figure 10. Model self-validation at natural cooling mode with various rotation speeds (a) Rear bearing (T1) temperature validation (b) Midpoint of housing (T2) temperature validation (c) Front bearing (T3) temperature validation (d) Midpoint of shaft (T4) temperature validation.
Sensors 18 00656 g010
Figure 11. External validation between measured temperature and predicted temperature of estimated TNM under stepwise rotational speed (3001, 5018, 7028 rpm).
Figure 11. External validation between measured temperature and predicted temperature of estimated TNM under stepwise rotational speed (3001, 5018, 7028 rpm).
Sensors 18 00656 g011
Figure 12. (a) Hankel singular values; (b) Short circuit time constant method with C2–C4 are short circuited; (c) Pole location of estimated TNM when comparing with the Model Order Reduction (MOR) model, SCTC model and varying with rotational speed.
Figure 12. (a) Hankel singular values; (b) Short circuit time constant method with C2–C4 are short circuited; (c) Pole location of estimated TNM when comparing with the Model Order Reduction (MOR) model, SCTC model and varying with rotational speed.
Sensors 18 00656 g012
Figure 13. Bode diagram of estimated TNM and 1st-order truncated model based on MOR.
Figure 13. Bode diagram of estimated TNM and 1st-order truncated model based on MOR.
Sensors 18 00656 g013
Figure 14. Schematic summary.
Figure 14. Schematic summary.
Sensors 18 00656 g014
Figure 15. Demonstration of the thermo-feature identification system.
Figure 15. Demonstration of the thermo-feature identification system.
Sensors 18 00656 g015
Table 1. BTSM Specification.
Table 1. BTSM Specification.
Sensor ElementAccuracy [°C]Resolution [°C]Measurement Range [°C]Power [mW]Module Size [mm3]
RTD±(0.1 + 0.0029|ϑ|)0.00489−40~1507Ø40×
Table 2. Estimated thermal parameters of TNM on operating mode.
Table 2. Estimated thermal parameters of TNM on operating mode.
ParameterValueParameterValueParameterValue
R12 [KW−1]0.1263Ra1 [KW−1]0.884Rav3 [KW−1]1.0032
R14 [KW−1]0.1537Ra2 [KW−1]1.393qf1 [W]37.162
R23 [KW−1]0.5954Ra3 [KW−1]3.071qf2 [W]28.406
R24 [KW−1]3.609Rav1 [KW−1]1.498qf3 [W]0.0043
R34 [KW−1]0.48Rav2 [KW−1]0.303qf4 [W]0.0073
Table 3. Estimated thermal parameters of TNM on natural cooling mode.
Table 3. Estimated thermal parameters of TNM on natural cooling mode.
ParameterValueParameterValueParameterValue
R12 [KW−1]0.1263Ra1 [KW−1]18.765C1 [JK−1]5375.8
R14 [KW−1]0.1537Ra2 [KW−1]1.393C2 [JK−1]3545.8
R23 [KW−1]0.5954Ra3 [KW−1]6.496C3 [JK−1]10,931.7
R24 [KW−1]6.737Ra4 [KW−1]2.497C4 [JK−1]625.4
R34 [KW−1]0.48

Share and Cite

MDPI and ACS Style

Lo, Y.-C.; Hu, Y.-C.; Chang, P.-Z. Parameter Estimation of the Thermal Network Model of a Machine Tool Spindle by Self-made Bluetooth Temperature Sensor Module. Sensors 2018, 18, 656. https://doi.org/10.3390/s18020656

AMA Style

Lo Y-C, Hu Y-C, Chang P-Z. Parameter Estimation of the Thermal Network Model of a Machine Tool Spindle by Self-made Bluetooth Temperature Sensor Module. Sensors. 2018; 18(2):656. https://doi.org/10.3390/s18020656

Chicago/Turabian Style

Lo, Yuan-Chieh, Yuh-Chung Hu, and Pei-Zen Chang. 2018. "Parameter Estimation of the Thermal Network Model of a Machine Tool Spindle by Self-made Bluetooth Temperature Sensor Module" Sensors 18, no. 2: 656. https://doi.org/10.3390/s18020656

APA Style

Lo, Y. -C., Hu, Y. -C., & Chang, P. -Z. (2018). Parameter Estimation of the Thermal Network Model of a Machine Tool Spindle by Self-made Bluetooth Temperature Sensor Module. Sensors, 18(2), 656. https://doi.org/10.3390/s18020656

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop