Abstract
The stable operation of a high-speed rotating rotor-bearing system is dependent on the internal damping of its materials. In this study, the dynamic behaviours of a rotor-shaft system with internal damping composite materials under the action of a temperature field are analysed. The temperature field will increase the tangential force generated by the internal damping of the composite material. The tangential force will also increase with the rotor speed, which can destabilise the rotor-shaft system. To better understand the dynamic behaviours of the system, we introduced a finite element calculation model of a rotor-shaft system based on a 3D high-order element (Solid186) to study the turbocharger rotor-bearing system in a temperature field. The analysis was done according to the modal damping coefficient, stability limit speed, and unbalance response. The results show that accurate prediction of internal damping energy dissipation in a temperature field is crucial for accurate prediction of rotor dynamic performance. This is an important step to understand dynamic rotor stress and rotor dynamic design.
1. Introduction
Turbochargers are mechanical devices that can improve fuel efficiency and reduce greenhouse gas emissions. The core component of a turbocharger is a rotor composed mainly of a turbine and a compressor. The turbine is crucial because it can recover energy from the exhaust gas and increase the intake air volume by driving the compressor [1]. It is characterized by light weight and high speed. To meet requirements for higher speeds, greater power density, and the ability to adapt to harsh operating environments for a long time, structural designers use various materials to manufacture the rotor parts. This introduces high requirements for the stable rotation of the rotor. Therefore, to analyse the stability of a rotor and its dynamic characteristics, it is necessary to accurately predict the damping effect. Damping is divided into external damping, such as bearing damping, and internal damping, such as material damping, which in turn is mainly modelled by viscous damping and hysteretic damping [2]. In a dynamic composite rotor, internal damping is meaningful only when the matrix can damp [3]. Also, most single materials, such as metals, have vibration damping characteristics similar to those of hysteretic internal damping but not viscous internal damping [4].
Many researchers have studied the effect of material damping on rotor dynamics and stability behaviour. Sin [5] studied the instability phenomenon of a composite shaft with internal damping and calculated the natural frequency and instability threshold of the shaft by using the finite element model for beams. They found that the stacking sequence of layers of composite materials, the arrangement directions of material fibers, and the transverse shear forces of materials affected the natural frequency and instability threshold of the shaft. Wittgren [6] studied the stability of flexible rotors with flexible supports. They found that the correct combination of the asymmetry of the internal damping of the shaft material and the anisotropy of bearings could significantly reduce the amplitude at critical speed. Montagnier [7] used Euler–Bernoulli shaft model to analyse the critical speed of a flexible internal resistance rotor supported by elastic bearings to obtain the stable working range of the rotor. They found that most of the rotor instability was caused by the internal damping of the material. Increasing the bearing damping could make the rotor run more stably, and they established a rotor instability analysis model that included the internal resistance of the material. Vitta [8] used linear and nonlinear evaluation methods to analyse the effect of internal damping on the dynamic behaviour of a shaft. The results show that the critical, subcritical, and stable limit speeds of a rotor can be obtained by linear evaluation. Newkirk [9] studied internal damping and found that a rotor may vibrate violently at a speed higher than the first critical speed. Genta [10] explained that the hysteretic damping of rotating structures is stable in the subcritical range but unstable in the supercritical range. All those researchers proved that the instability speed is higher than the first-order critical speed when considering the hysteretic damping of rotating structures. Also, many researchers have studied the combined effects of material damping and bearing damping. The results show that increasing the bearing damping can improve the stability of the rotor, but increasing the material damping can reduce the instability threshold [6].
The effect of the surrounding temperature on material damping of a turbocharger rotor cannot be disregarded in high-temperature environments. It is very important to understand the dynamic behaviours of the structure in various temperatures when designing a rotor to operate in thermal extremes. San [11] established a nonlinear rotor-bearing finite element model. The model considered the thermal effects of lubricating oil and the thermal expansions of the rotor shaft and bearing, and San Andres analysed the natural frequency, stability, and unbalanced response of the rotor. The results show that the natural frequency and synchronous speed amplitude obtained by the simulation were completely consistent with the experimental values, which verified the feasibility of the finite element model. Zych [12] used a finite element method to calculate the thermal stress of a radial axial microturbine in a high-temperature environment. In the calculation, they considered the mass of the disc, the rotation speed of the rotor, and the complex shape at the rear of the disc. Their work showed that numerical calculation helps to choose the best optimization method, and they reduced the turbine’s von Mises stress by approximately 45%. Jeyaraj [13] uniformly heated a thin plate that incorporated various internal damping composite materials and did free vibration and forced vibration analyses. The study found that, with increasing temperature, the vibration amplitude (response) of the thin plate structure decreased, but the modal loss factor increased significantly with increasing temperature, thereby reducing the vibration amplitude. Its response frequency also decreased with increasing temperature. Guo [14] did thermal vibration analysis of a rotating beam structure with constrained layer damping. The study found that, as the temperature increased, the modal frequency of the beam structure decreased accordingly, and the damping ratio increased accordingly. That study provides a basis for dynamic analysis of high-speed rotating blades in various thermal environments. However, the number of articles on this subject is limited.
Therefore, this research studied the dynamic characteristics of internal material damping and an oil film force turbocharger rotor-bearing system under thermal environments (temperature fields). We used the conjugate heat transfer (CHT) method to simulate the temperature field of the solid part of a rotor-shaft system. The temperature field was coupled with the finite element model of the rotor, compared with the instance without considering the temperature field. The rotor finite element is verified by experiment.
2. Numerical Model and Research Methods
2.1. Fluid and Heat Transfer Sections
Before the thermal modal analysis, the aerodynamic thermal analysis was performed. The current industry standard modelling approaches assume the turbine and compressor operate under adiabatic conditions [15]. The CHT simulation method [16–22] is used to obtain the temperature distribution of the rotor, and the node temperature is provided for the modal simulation part.
2.1.1. Turbocharge and Compressor Geometry
The turbocharger had an impeller with 10 blades, and the compressor had an impeller with 6 blades and 6 splitters. The turbocharger and compressor parameters are shown in Table 1. To improve the accuracy of the calculations, we modelled both the turbocharger and the compressor using full passages. The solid parts and the air flow are shown in Figure 1.
2.1.2. Grid Generation
CHT involves the direct coupling of fluids and solids. ICEM grid discretisation software (Ansys, Inc., USA) uses the same numerical principles and grid discretisation for both regions. This allows the noninterpolated exchange of heat flux between adjacent grids [23]. The calculation accuracy of CHT is very sensitive to the resolution of the fluid boundary layer grid. Therefore, the dimensionless distance or less (in this paper ) of the wall distance of the first layer of the grid can determine the local heat flux with enough accuracy. As shown in Figures 2 and 3, the fluid domain, solid domain (rotor impeller), and boundary layer grid were generated for CHT calculations. The total number of global grids was 12,329,631, with 2,067,903 global grid nodes. Following the shape of the outer contour of the rotor, we used a triangular tetrahedral mesh with good adaptability to the outer contour, and we used smoothing to optimize the mesh. The impeller grid and the fluid grid were connected in a general grid interface mode in Ansys CFX software.
2.1.3. Boundary Conditions
We determined the boundary conditions for aerodynamic thermal analysis using experimental data provided by the turbocharger company; the boundary conditions are shown in Table 2.
The heat transfer of the surrounding air was disregarded, and the out-wall of the turbocharger was assumed to be adiabatic. We applied “no-slip” boundary conditions to all inner walls. An interface was added between the rotating domain and the fixed domain, and the interface was connected by a “frozen rotor” [24].
2.1.4. Numerical Methodology
CHT refers to a coupled heat transfer phenomenon in which the thermal properties of two materials occur through a medium or in direct contact. The CHT method can calculate the heat transfer between fluid and solid and calculate the temperatures of fluids and solids at the same time. In this study, we used the commercial software Ansys CFX for numerical simulation. CFX is a computational fluid dynamics software package based on the control volume method to solve Navier–Stokes equations. In the fluid domain, the mass conservation, momentum, and energy transport equations are described aswhere , , , , , , , and represent the velocity vector, density, pressure, stress tensor, thermal conductivity, temperature, total enthalpy, and time of the fluid, respectively.
In a solid domain, the conservation of energy equation can explain the heat transfer caused by solid motion, conduction, and a volume heat source. The energy equation iswhere is the velocity vector of the solid, ; is the optional volume heat source, ; is the density of the solid; and is the thermal conductivity of the solid.
This study did not directly solve (1)–(4). Instead, they were converted to the steady-state Reynolds average Navier–Stokes method to calculate turbulence. The fluid medium (exhaust gas and air) is an ideal gas, and the shear stress transmission (SST) turbulence model was used because the SST turbulence model has good accuracy for CHT calculations [20]. That model has both the accuracy of the model in high-pressure gradient flow boundary layer prediction and the stability of the model in mainstream prediction [18].
2.2. Modal Part
2.2.1. Rotor-Shaft System Geometry
The model of the rotor-shaft system was provided by the turbocharger company and agreed with the real one. The rotor-shaft system comprised a turbine wheel, a compressor wheel, a rotating shaft, a floating ring bearing, a thrust bearing, a seal sleeve, and a nut, as Figure 4 shows. The turbine wheel and shaft were connected by friction welding, whereas the compressor was axially and circumferentially fixed by a left-hand nut.
The materials of the main parts of the rotor-shaft system are shown in Table 3, based on the actual situation.
This study investigated the effect of temperature and found that the material properties changed as the temperature changed. The function curves of the material properties of each part as a function of temperature are shown in Figure 5 [25, 26].
(a)
(b)
(c)
(d)
Publication [29] proposes a set of damping coefficients, and we verified the influence of damping coefficients on the stability limit displacement of the rotor through theoretical analysis and numerical simulation. We selected this set of damping coefficients based on the theory of the publication.
Table 4 shows what we assumed was the material damping coefficient of each part of the rotor.
2.2.2. Floating Ring Bearing Parameters
The rotor shaft was supported by a floating ring bearing whose detailed parameters are shown in Table 5.
To aid the simulation, we combined the stiffness and damping coefficients of the inner and outer oil film into a total impedance (equivalent stiffness and damping coefficients) according to [27, 28]. Combining the stiffness and damping coefficients of the inner and outer oil film provided by the turbocharger company with the parameters of the floating ring bearing, we calculated the total impedance using Dyrobes rotor dynamics software (Dyrobes, USA). Figure 6 shows the total impedance (equivalent stiffness and damping coefficients) of the floating ring bearing at various speeds.
(a)
(b)
2.2.3. Finite Element Model of Rotor
This study used a large-scale general software Ansys Workbench for simulation. To obtain the motion trace on the axis, the rotor shaft was equally divided into four parts along the axis direction, as shown in Figure 7. The blue lines (solid and dotted) in the figure represent the division positions.
The solid model was discretised, and an Ansys Solid186 element was used. The Solid186 element is a high-order, 3D 20-node element with quadratic displacement characteristics. The element is defined by 20 nodes, each node with three degrees of freedom (translation of each node in the x-, y-, and z-directions). The element supports plasticity, superelasticity, creep, stress hardening, large deflection, and large strain capacity. In the shaft part, the grid was generated by a sweep, as shown in Figure 8, depicting the “” element after segmentation.
We transformed (deformed) the original hexahedron element to an approximately one-quarter cylinder element. In Figure 8, the red dotted line represents the actual edge of the element, and the red solid line represents the edge of the element. The advantage of a solid model is that it retains the properties of the original element and describes the outer contour of the shaft model well.
To adapt to the more complicated outer contour of the blade, tetrahedral, pyramid, and prism methods were used to generate the unit in the blade and the remaining part. Figure 9 shows the (), (), and () elements generated by three methods.
The effects of rotational inertia, translational inertia, gyro moment, support stiffness, support damping, material damping, rotational damping, and thermal stress stiffness were considered in the model to establish the equation of motion. The rotor model was discretised into elements, and the total number of nodes became . The node displacement vector of the element iswhere is the node number and the superscript is the matrix/vector transpose. These writing formats are also mentioned later in the text.
After combining the governing equations of all elements and combining the boundary conditions, the equation of motion of the rotor-bearing system iswhere ; ; and .
Here, is a mass matrix that includes mass components caused by translation along the axis () and rotation along the axis (). is a damping matrix (asymmetric matrix coefficient of the velocity vector) that includes a skew-symmetric gyro matrix (), a constant structural damping matrix () with hysteresis characteristics caused by internal friction of the material and damping caused by the floating ring bearing support matrix (). is the stiffness matrix (matrix coefficient of the displacement vector) that includes the stiffness matrix () caused by the support of the floating ring bearing, the spin-softening bending matrix () generated by the rotor rotation, the stiffness matrix () of the elastic modulus, the stiffness matrix () caused by the nonuniform thermal stress, and thermal deformation caused by temperature changes. The size of all matrices is . is the rotor speed in rpm. is the material damping (Rayleigh damping), stiffness damping coefficient.
The total displacement and global force vectors and are expressed as
Writing (6) in state-space form yieldswhere matrices and ; state vector ; and force vector arewhere is an empty matrix of size and an empty vector of size .
Figure 10 shows a simplified finite element model of the rotor-shaft system. The grey-blue points in the figure are used to observe the modal shapes, line trajectories, and unbalanced harmonic response analysis.
Figure 11 shows the rotor finite element global grid. The total number of grids was 3,176,485, with 546,067 grid nodes. In the geometry software, each portion of the rotor was treated as a part to avoid extra contact stiffness in the finite element simulation and prevent affecting the simulation accuracy.
To aid the comparison, we considered two cases in this study: case 1, the rotor-shaft system without the influence of temperature, and case 2, the rotor-shaft system affected by temperature. In the following sections, r-th forward and r-th backward modes are represented as “r-F” and “r-B” modes, respectively.
3. Results and Discussion
3.1. Heat Transfer Results and Experimental Verification
After the calculation, we verified the numerical results of the mass flow and the supercharger ratio (expansion ratio). The numerical results agreed with the experimental data, and the error control was approximately 7%, as Figure 12 shows.
(a)
(b)
The temperature distribution of the rotor was obtained through postprocessing, as Figure 13 shows; the temperature distribution law was the same as that in the Bohn [19] study result and is shown in Figure 14.
Figures 13 and 14 show that the highest temperature was at the outer edge of the turbine wheel blades and that the lowest temperature was at the compressor wheel. Overall, the rotor had a large temperature gradient.
3.2. Modal Results
3.2.1. Campbell Diagram and Stability Diagram of Normal Temperature Environment
Figure 15 shows the Campbell diagram of the rotor-shaft system in case 1, that is, without considering the temperature. The figure was obtained from the whirl frequency (obtained from the imaginary part of the eigenvalue) and the rotational speed. To present a clear figure and express its basic characteristics, only the first four whirling directions are described. In the figure, there are two forward whirls (the whirl direction of the rotor agrees with the rotation direction of the speed), which are marked as 1F and 2F in order of increased frequency. There are also two backward whirls (the whirl direction of the rotor is opposite to the rotation direction of the speed) that are marked as 3B and 4B, respectively. The synchronous whirl line is called the synchronous excitation line. In the figure, the velocities corresponding to the intersection of the synchronous whirl line and the frequency line of the forward whirl are marked as and , respectively, and were calculated as follows: and . Because the backward whirl could not be excited by the unbalanced force, it is not marked or considered.
When the speed was in the critical speed range between and in the 1F mode, the whirling frequency had a clear increasing trend. That trend may have been caused by floating ring bearings, because in the critical speed range of and , the stiffness and damping of the oil film changed sharply. When the critical speed of was exceeded, the stiffness and damping of the oil film changed smoothly, and the eddy frequency in the corresponding 1F mode also increased at a constant rate.
Figure 16 shows that the curve of the damping ratio of the four modes in case 1 changed with the rotation speed. The modal damping ratio is the ratio of actual damping to the critical damping of a certain mode, and the damping ratio is the decisive parameter to describe the stability. A positive damping ratio means that the rotor-shaft system was in a stable region when the vibration energy was dissipated. A negative damping ratio means that the rotation energy supported the rotor rotation by adding vibration energy, resulting in the rotor-shaft system becoming unstable. The point where the damping ratio curve intersects with the zero line is the stable limit speed (SLS). The figure shows that the rotor-shaft system became unstable when in 1F mode and then became unstable when in 2F mode.
Figure 17 shows the four modes close to the SLS 1 stability limit speed at 55,000 rpm. The black dotted line represents the axis, the combination of the orange line and the orange sphere represents the starting position of the trajectory, and the incomplete track curve represents the whirl direction. The figure shows that the 1F mode was cylindrical, thus combining rigid body motion with rotor bending. Furthermore, the rotor was close to the unstable state at this time. In contrast, the 2F mode was a conical, bending mode. The figure shows that the 3B mode was conical, with both rigid body motion and rotor bending. Lastly, the 4B mode was cylindrical, with both rigid body motion and rotor bending.
(a)
(b)
(c)
(d)
Figure 18 shows the four modes close to the SLS 2 stability limit speed at 97,946 rpm. The rotor was close to the unstable state at this time. The 1F and 3B modes were cylindrical and conical, combining rigid body motion and rotor bending. In contrast, the 2F mode was a conical, bending mode.
(a)
(b)
(c)
(d)
3.2.2. Effect of Temperature on the Campbell Diagram and Stability Diagram
We exported the temperature data of each node in Figure 12 and inserted the temperature value into the corresponding rotor-shaft system node through the Ansys software commands function. We analysed the critical speed and stability of the rotor-shaft system in the thermal environment.
Figure 19 shows that temperature affected the intrinsic frequency (whirl frequency). Compared with case 1, the whirl-frequency curves in the 1F, 2F, and 3B modes were slightly lower, whereas for the 4B mode, the whirl-frequency curve was much lower. In the figure, the starting point at the left end (when the rotating speed was 0 rpm) was taken as the observation object. In case 1, the frequency was 2,300.9 Hz, whereas in case 2, the frequency decreased to 2,182.5 Hz. Because the temperature changed the material’s elastic modulus and Poisson’s ratio, it affected its structural stiffness. For higher-order frequencies, the effect on frequency became greater as the order increased. However, the effect on the critical speed, , did not change much, whereas decreased from 43,622 rpm to 41,815 rpm. This decrease was also caused by the intrinsic frequency drop.
Figure 20 shows the stability diagram in case 2 and that the temperature had an important effect on the damping ratio and SLS. Compared with case 1, SLS 1 of the 1F mode decreased from the original 56,550 rpm to 55,078 rpm and made the 1F mode enter the unstable region ahead of time. The 2F mode change was more obvious, and the same SLS 2 decreased from the original 97,946 rpm to 95,969 rpm. The reason for this can be explained by the following formula: . Here, is the damping ratio coefficient of the mode, is the mass damping coefficient, is the stiffness damping coefficient, and is the natural frequency of the mode. In most cases, () can be disregarded; that is, . did not change, whereas decreased because of the effect of temperature, resulting in a decrease in . In other words, SLS 1 and SLS 2 in the 1F and 2F modes were reduced. Therefore, when predicting the stability limit reset and dynamic performance of a turbocharger rotor, the effect of temperature should have a great influence on the dynamic design of the rotor.
3.2.3. Unbalance Response Analysis
We applied a dummy unbalance mass of 0.01 at node 2, and the frequency range of the steady-state synchronous response analysis was (corresponding speed range ). We then measured the response amplitude of the turbine end floating ring bearing center position node 5. In Figure 21, case 1 is represented by a black curve (normal temperature), and case 2 by a red curve (high temperature). In case 1 (case 2), the curves represent 6F, 6B, and 20F modes. The frequency at the wave crest corresponding to the frequency corresponding to the fast reverse rotation (the frequency of the rotation speed within one minute) in the case 1 (case 2) Campbell diagram, as shown in Figures 22 and 23), is similar to the black (red) number in Figure 20. Because the steering mode is dominant in the response diagram, the critical speeds and in case 1 (case 2) 1F and 2F modes are not shown in the diagram and followed a phenomenon similar to that reported by Chouksey and Roy [29, 30]. Because of the effect of temperature, the steering frequency of case 2 shifted to the left, making the change of response amplitude more obvious. The response value at the peak of the 6F and 6B modes decreased from 0.147 mm to 0.014 mm, and that at the peak of the 20F and 20B modes decreased from 0.036 mm to 0.0037 mm.
3.2.4. Verification of Rotor FEM
Natural frequency analysis is very important as the basis of a dynamic analysis. We obtained the first two natural frequencies of the rotor by hammering experiments and verified them by numerical simulation. Figure 24 shows the measured parts (rotor system), force hammer, and piezoelectric acceleration sensor. To reduce the measurement error rate, we adopted a sensor with a nonlinearity of less than 1.2% and a sensitivity of 10 mv/g and placed the sensor in the middle of the rotor shaft. Table 6 compares the natural frequency experimental value with the simulated value. The experimental and simulated values of the first natural frequency were 1386 Hz and 1453 Hz, respectively. The second-order natural frequencies were 4100 Hz and 4383 Hz, respectively, and the error rate was within 7%. But, the higher mass of the accelerometer is the main cause of error value between the experimental value and the simulated value.
A more favorable verification was the oil film data shown in Figure 6. The oil film data were experimental data provided by the turbocharger company. In the simulation, the oil film was the key factor in supporting or restraining the rotor. This also illustrates the reliability of the simulation method.
4. Conclusions
This study used the CHT numerical simulation method to obtain the temperature of a rotor-shaft system and found that there was a large temperature gradient when the system was in use.
Splitting the rotation shaft solved the problem that a solid model cannot obtain an axis orbit on the Ansys Workbench software platform.
Analysing the Campbell diagram shows that the temperature reduced the stiffness of the rotor structure, causing the intrinsic frequency to decrease, especially for higher-order frequencies. The decrease in stiffness was greater and the critical speed (Ωcr1) in the 1F mode did not change much, whereas the critical speed (Ωcr2) in the 2F mode decreased a large amount. Therefore, dynamic design and stability prediction of a turbocharger rotor must consider the temperature during its operation.
For the analysis of an unbalanced response, the response corresponding to the critical speed in the 1F and 2F modes does not appear in Figure 21 because these two modes have higher damping and because the two steering modes were dominant. However, the response diagram can well reflect the rotor steering mode. Also, the temperature effect on the response was mainly reflected in the response amplitude, and the response frequency was slightly shifted to the left (reduced), which was similar to the decreased intrinsic frequency relation.
Abbreviations
| CHT: | Conjugate heat transfer |
| SWL: | Synchronous whirl line |
| SLS: | Stable limit speed |
| TC: | Turbocharger |
| r-F: | r-th forward whirl mode |
| r-B: | r-th backward whirl mode. |
| : | Mass damping coefficient |
| : | Stiffness damping coefficient |
| : | Natural base |
| : | Stress |
| : | Damping ratio coefficient |
| : | Dimensionless wall thickness. |
| brg: | Bearing |
| cr: | Critical |
| con: | Constant |
| e: | Element |
| E: | Energy |
| f: | Fluid |
| h: | Harmonic |
| s: | Solid |
| sh: | Shaft |
| trs: | Translation |
| rot: | Rotation |
| tem: | Temperature |
| tot: | Total. |
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This work was supported by the Key Research & Developmental Program of Shandong Province (Grant no. 2019GGX104104). The authors are grateful for technical support and the experimental data from Kangyue Technology Co., Ltd.