Coupled Fluid-Thermal Analysis for Induction Motors with Broken Bars Operating under the Rated Load

Abstract: Thermal stress of the rotor in a squirrel cage induction motor is generated due to the temperature rise, it is also one of the factors causing the broken bar fault because the structure of the rotor would be destroyed if the stress of the rotor bars exceed the strength limit. The coupled fluid-thermal analysis for the induction motor with healthy and broken bar rotors is performed in this paper. Much concern has been committed to establishment of the fluid model on the basis of computational fluid dynamic (CFD) theory. The heat field of the prototypes is analysed so that the effect of the asymmetrical rotor on the motor heat performance can be investigated in depth. Eventually, the efficiency of the presented model and method, for the totally enclosed fan cooled (TEFC) induction motor, can be verified through experimental results. In addition, this paper reports a quantitative analysis of the heat flux distribution of the fault rotor, and the heat flux density of the bars is investigated in detail. Then, the part most likely to break in the rotor as a result of the thermal load is identified.


Introduction
As a key performance parameter of electrical machines, the temperature rise has been of wide concern in the industrial field.The temperature rise can limit the rated power and reduce the motor efficiency, and it also imposes special demands on the material of the motor.Therefore, in the process of the motor design, the temperature limits should be calculated and checked for satisfying the requirements of functionality and safe.If the operating temperature exceeds the insulation material working limit of the motor [1], it will cause an accelerated aging effect or even failure of insulation [2,3], and then machine life expectance can be seriously affected [4].Accordingly, it is vital to plot the thermal map of the electric machine, working in the most severe condition.
As well known, the accuracy of thermal analysis greatly depends on the heat transfer coefficient on heat transfer surfaces, so it is essential for designers to get an accurate analysis of the fluid field.However, inner and outer environment of the motor is complex, in some conventional methods [5][6][7], the fluid model is usually simplified, or the fluid velocity on heat transfer surfaces is regarded as constant, so the corresponding heat transfer coefficient is also constant.These studies provide good ways to fluid field analysis, but the continuity of the fluid field distribution is not included.In view of this problem, reference [8] proposed a coupled calculation of fluid and thermal field for high speed motor considering tooth slots effects, and [9] proposed analysis of the electromagnetic-fluid-thermal field in a permanent magnet linear motor.The permanent magnet motor with high speed of 60,000 r/min is discussed in [10], according to computational fluid dynamic (CFD) theory, the parameters of the equivalent thermal circuit can be obtained.For an air-cooled disc type electrical machine, the measurement and simulation of stator heat transfer were performed in [11] and the CFD analysis was introduced.The CFD (computational fluid dynamic) analysis and the thermal equivalent circuit model have been established for a totally enclosed fan cooled (TEFC) induction motor, and the simulation data are compared with experimental data measured by sensors [12].Moreover, many fluid analyses of the motor are based on local fluid field, very little research has addressed the complete fluid space, especially for small and medium size induction motors.
In addition, apparent thermal stresses can be observed with the temperature rise of motors.It frequently leads to cage fracture in the joint of the bars and the rings and accompanying asymmetric rotor operation, which is a serious accident in the running state of the motor.The thermal field distribution research on asymmetry rotor can provide some references for the prevention and diagnosis of fault.With the increase of rotor bar failures, it has received some attention from researchers.In [13], an innovative method of fault diagnosis for the induction motor with broken bars fault is proposed.A new methodology was introduced in [14], which integrates vibration and current analysis by associating signal spectrum analysis to improve the detectability for the motor, under the mechanically loaded and unloaded running state respectively.Reference [15] evaluated the internal faults of the induction motor with the squirrel cage by the finite element method as well as discrete wavelet transform.Reference [16] presented an experimental research of the fault rotor motor, regarding time frequency evolution effects, and the fault factor can be obtained.The detection methods for induction motor with rotor bar fault are proposed in [17,18], and different diagnostic mediums are compared in two studies.The purpose of [19] is to report the effect of partially fault rotor bar based on experimental research in induction motor.The running characteristics of the motor under different operating conditions are obtained by discrete wavelet transform analysis.The effect of saturation level due to rotor bars fault in induction motor was introduced in [20], by means of fundamental electromagnetic laws and finite element modelling.Obviously, the existing studies works relating to thermal analysis of the motor with broken bar is still much less than those associated with other aspects.
Considering the aforementioned issues, in this paper, the coupled fluid-thermal analysis of the prototype running under normal and faulty conditions is performed.To effectively simulate the fluid field, the 3D model for the prototype motor and the fluid domains are built, and the coupled model is solved by CFD simulation.The software "ANSYS" was used to analyze the induction motor with healthy and faulty rotor.The 3D geometry model construction of the induction motor is studied by commercial software ANSYS ICEM (15, ANSYS Inc., Canonsburg, PA, USA) programs, and the calculation of the coupled fluid-thermal field, including meshing, solving and post-processing is achieved by commercial software ANSYS FLUENT (15, ANSYS Inc., Canonsburg, PA, USA).In the stage of pre-processing, the outer approximate infinite fluid field and 3D finite element construction model are established, with the ventilation structure and characteristics of heat transfer taken into account.In the light of Newton's heat transfer law and Fourier's theorem, the heat transfer coefficients of the motor can be obtained, and the steady thermal field of the motor running at full capacity is computed.The calculated results are close to the experimental data, and the analytical method in the paper proved to be valuable can provide references for the other middle and small induction motors.

Basic Parameters of the Motor and Laboratory Testing Platform
The major data of the prototype motor are listed in Table 1.The air velocity of fluid field and temperature of the motor are also tested in the laboratory under healthy and broken bar conditions.The test flat has been shown in Figure 1, and the same rotors, including one healthy rotor and two faulty rotor are prepared for the test purpose.The motor used for experimental study is a three-phase induction motor, 50 Hz, 2 poles, 1.1 kW, and the healthy squired cage rotor is made of 16 bars.The air velocity along the axial direction on shell surface can be accessed by a hot wire anemometer during the experiment as shown in Figure 1.In order to evaluate the heat behaviour of the induction motor, a complete thermal evaluation of the prototype motor fitted with thermistor is carried out in the laboratory, so the temperature in all the spots of interest can be obtained.The temperature on the motor shell is measured by an infrared thermometry probe.

CFD Model and Temperature Measured Positions
The 3D model of the prototype motor is employed in the paper, as shown in Figure 2. Meanwhile the presented model is intended to calculate the fluid and temperature field.The structured grids and unstructured grid are adopted in the paper.The motor is meshed by structured grid due to complex structures, and the external fluid space is meshed by the unstructured grids.There are 3,160,630 mesh cells in the simulations.The air velocity along the axial direction on shell surface can be accessed by a hot wire anemometer during the experiment as shown in Figure 1.In order to evaluate the heat behaviour of the induction motor, a complete thermal evaluation of the prototype motor fitted with thermistor is carried out in the laboratory, so the temperature in all the spots of interest can be obtained.The temperature on the motor shell is measured by an infrared thermometry probe.

CFD Model and Temperature Measured Positions
The 3D model of the prototype motor is employed in the paper, as shown in Figure 2. Meanwhile the presented model is intended to calculate the fluid and temperature field.The structured grids and unstructured grid are adopted in the paper.The motor is meshed by structured grid due to complex structures, and the external fluid space is meshed by the unstructured grids.There are 3,160,630 mesh cells in the simulations.The air velocity along the axial direction on shell surface can be accessed by a hot wire anemometer during the experiment as shown in Figure 1.In order to evaluate the heat behaviour of the induction motor, a complete thermal evaluation of the prototype motor fitted with thermistor is carried out in the laboratory, so the temperature in all the spots of interest can be obtained.The temperature on the motor shell is measured by an infrared thermometry probe.

CFD Model and Temperature Measured Positions
The 3D model of the prototype motor is employed in the paper, as shown in Figure 2. Meanwhile the presented model is intended to calculate the fluid and temperature field.The structured grids and unstructured grid are adopted in the paper.The motor is meshed by structured grid due to complex structures, and the external fluid space is meshed by the unstructured grids.There are 3,160,630 mesh cells in the simulations.Due to the irregular arrangement of the winding in the stator slots, the model must be simplified based on the following hypotheses [21]: 1.
The stator winding has a good dipping process.

2.
Insulated lacquer on the stator winding surface is uniformly distributed.

3.
The equivalent insulation layer in the stator slots is tightly connected with the stator core.4.
The thermal difference of all winding in a stator slot can be neglected.
Based on the above hypothesis, all copper wires can be equivalent treated as a copper bar placed in the stator slot center.The insulating paper, impregnating varnish, varnish film of the copper and the air between wires are dealt with an equivalent insulating layer.The equivalent method is shown in the Figures 3 and 4. The equivalent heat conduction coefficient of the insulation layer can be calculated by Equation (1) [21].
where, λ eq is equivalent heat conduction coefficient of the insulation layer in the stator slots, δ i is equivalent thickness of each insulating material (δ 1 is thickness of the slot wedge, the value is 0.5 mm, δ 2 is thickness of the air gap between slot wedge and stator core, the value is 0.15 mm, δ 3 is thickness of the enameled wire insulation and impregnating varnish, the value is 0.96 mm), λ i is average thermal conductivity of each material.
Energies 2018, 11, x FOR PEER REVIEW 4 of 17 Due to the irregular arrangement of the winding in the stator slots, the model must be simplified based on the following hypotheses [21]: 1.The stator winding has a good dipping process.2. Insulated lacquer on the stator winding surface is uniformly distributed.3. The equivalent insulation layer in the stator slots is tightly connected with the stator core.4. The thermal difference of all winding in a stator slot can be neglected.
Based on the above hypothesis, all copper wires can be equivalent treated as a copper bar placed in the stator slot center.The insulating paper, impregnating varnish, varnish film of the copper and the air between wires are dealt with an equivalent insulating layer.The equivalent method is shown in the Figures 3 and 4. The equivalent heat conduction coefficient of the insulation layer can be calculated by Equation (1) [21].
( ) where, λeq is equivalent heat conduction coefficient of the insulation layer in the stator slots, δi is equivalent thickness of each insulating material (δ1 is thickness of the slot wedge, the value is 0.5 mm, δ2 is thickness of the air gap between slot wedge and stator core, the value is 0.15 mm, δ3 is thickness of the enameled wire insulation and impregnating varnish, the value is 0.96 mm), λi is average thermal conductivity of each material.Furthermore, because of the low heat conduction coefficient of air, the thermal contact resistance of the assembly space between the stator and the frame cannot negligible, and it has a great effect on the thermal analysis of the motor.The thermal contact resistance is defined as follows [22]: where, Rc is the thermal contact resistance, δ is the gap length between the stator and motor frame, the value is 0.5 × 10 −4 m, λ is the heat conduction coefficient of the air, Aa is the contact area between the stator and motor frame.The material properties are isotropic, and the list of material properties is given in Table 2. Due to the irregular arrangement of the winding in the stator slots, the model must be simplified based on the following hypotheses [21]: 1.The stator winding has a good dipping process.2. Insulated lacquer on the stator winding surface is uniformly distributed.3. The equivalent insulation layer in the stator slots is tightly connected with the stator core.4. The thermal difference of all winding in a stator slot can be neglected.
Based on the above hypothesis, all copper wires can be equivalent treated as a copper bar placed in the stator slot center.The insulating paper, impregnating varnish, varnish film of the copper and the air between wires are dealt with an equivalent insulating layer.The equivalent method is shown in the Figures 3 and 4. The equivalent heat conduction coefficient of the insulation layer can be calculated by Equation (1) [21].
( ) where, λeq is equivalent heat conduction coefficient of the insulation layer in the stator slots, δi is equivalent thickness of each insulating material (δ1 is thickness of the slot wedge, the value is 0.5 mm, δ2 is thickness of the air gap between slot wedge and stator core, the value is 0.15 mm, δ3 is thickness of the enameled wire insulation and impregnating varnish, the value is 0.96 mm), λi is average thermal conductivity of each material.Furthermore, because of the low heat conduction coefficient of air, the thermal contact resistance of the assembly space between the stator and the frame cannot negligible, and it has a great effect on the thermal analysis of the motor.The thermal contact resistance is defined as follows [22]: where, Rc is the thermal contact resistance, δ is the gap length between the stator and motor frame, the value is 0.5 × 10 −4 m, λ is the heat conduction coefficient of the air, Aa is the contact area between the stator and motor frame.The material properties are isotropic, and the list of material properties is given in Table 2. Furthermore, because of the low heat conduction coefficient of air, the thermal contact resistance of the assembly space between the stator and the frame cannot negligible, and it has a great effect on the thermal analysis of the motor.The thermal contact resistance is defined as follows [22]: where, R c is the thermal contact resistance, δ is the gap length between the stator and motor frame, the value is 0.5 × 10 −4 m, λ is the heat conduction coefficient of the air, A a is the contact area between the stator and motor frame.The material properties are isotropic, and the list of material properties is given in Table 2.The test points of the temperature A to G given in Figure 2, which generally report the temperature rise of the motor.The rotor temperature can be obtained from the sensor A planted in rotor bar.The aim of sensor B is to get the air-gap temperature, but actually it can be only mounted on stator wedge to obtain temperature of the air-gap near the stator.The sensors C and D are inserted inside two stator slots center respectively, and one is near the junction box and the other is placed in the opposite side of the junction box.E, F, G are three measured points of the motor frame.It should be noted that all test positions are located in the middle of the motor along the axial direction.

Fundamental Equations of the Fluid and Thermal Field
The fundamental equations governing fluid flow are solved for calculating the fluid field in this paper.As the range of the Reynolds number is greater than 2300, the turbulence equation is also required to describe the fluid flow.According to the CFD theory, the standard k-ε equation can be written as follows [23]: where, C µ , C 1ε and C 2ε are constants, Besides, 3D heat transfer equations are required for the temperature calculation [24,25]: where T is the temperature, T f is the fluid temperature, q is the heat generation per unit volume, α is the surface heat transfer coefficient, λ x , λ y , λ z are the coefficient of thermal conductivity in x, y, z direction respectively, n is the unit normal vector of the solution domain boundary, S j and S S are adiabatic surface and heat-dissipating surface respectively.The initial data α is constant, and it is variable on the heat transfer surface when numerical solutions are convergent.

Fundamental Assumptions
To simplify calculation, some reasonable assumptions are proposed as follows [26]: 1.
Thermal resistance of windings is neglected because of its good heat conducting properties.

2.
The thermal contact resistances between the rotor core and rotor bars is negligible.

3.
Stray loss is focused on the tooth tip of stator and rotor.

Boundary Conditions
The coupled fluid-thermal model are adopted in this paper, and the thermal boundary conditions of heat transfer surfaces can be got by the fluid field calculation.Compared to traditional empirical equations, the temperature results obtained through CFD simulation are more accurate.
In this section, the fluid domains of the motor were specially considered, so that they are closer to the real fact, as shown in Figure 5, the motor is placed in an approximate infinite fluid space.Considering the calculation time and accuracy, and the 20 times volume of the motor was chosen, it is appropriate for the computational domain.The external fluid space consists of two parts: one is the rotating fluid domain that is formed by a 360 degree sweep of the fan blade (see Figure 6) and the other is the rest of the fluid space.
adiabatic surface and heat-dissipating surface respectively.The initial data α is constant, and it is variable on the heat transfer surface when numerical solutions are convergent.

Fundamental Assumptions
To simplify calculation, some reasonable assumptions are proposed as follows [26]: 1. Thermal resistance of windings is neglected because of its good heat conducting properties.2. The thermal contact resistances between the rotor core and rotor bars is negligible.3. Stray loss is focused on the tooth tip of stator and rotor.

Boundary Conditions
The coupled fluid-thermal model are adopted in this paper, and the thermal boundary conditions of heat transfer surfaces can be got by the fluid field calculation.Compared to traditional empirical equations, the temperature results obtained through CFD simulation are more accurate.
In this section, the fluid domains of the motor were specially considered, so that they are closer to the real fact, as shown in Figure 5, the motor is placed in an approximate infinite fluid space.Considering the calculation time and accuracy, and the 20 times volume of the motor was chosen, it is appropriate for the computational domain.The external fluid space consists of two parts: one is the rotating fluid domain that is formed by a 360 degree sweep of the fan blade (see Figure 6) and the other is the rest of the fluid space.adiabatic surface and heat-dissipating surface respectively.The initial data α is constant, and it is variable on the heat transfer surface when numerical solutions are convergent.

Fundamental Assumptions
To simplify calculation, some reasonable assumptions are proposed as follows [26]: 1. Thermal resistance of windings is neglected because of its good heat conducting properties.
2. The thermal contact resistances between the rotor core and rotor bars is negligible.
3. Stray loss is focused on the tooth tip of stator and rotor.

Boundary Conditions
The coupled fluid-thermal model are adopted in this paper, and the thermal boundary conditions of heat transfer surfaces can be got by the fluid field calculation.Compared to traditional empirical equations, the temperature results obtained through CFD simulation are more accurate.
In this section, the fluid domains of the motor were specially considered, so that they are closer to the real fact, as shown in Figure 5, the motor is placed in an approximate infinite fluid space.Considering the calculation time and accuracy, and the 20 times volume of the motor was chosen, it is appropriate for the computational domain.The external fluid space consists of two parts: one is the rotating fluid domain that is formed by a 360 degree sweep of the fan blade (see Figure 6) and the other is the rest of the fluid space.Figure 7 is the section diagram of air-gap fluid field, and it is divided into two regions, including the rotary and stationary fluid region.The rotary fluid region is close to the rotor, and the stationary fluid region is close to the stator.On the basis of the reference coordinate system, the motive fluid domain and the rotor rotate synchronously, and the speed is given as a relative speed.The detailed boundary conditions of the fluid field are given in Table 3. motive fluid domain and the rotor rotate synchronously, and the speed is given as a relative speed.
The detailed boundary conditions of the fluid field are given in Table 3.

Losses
As well known the losses, in an induction motor bring about the temperature rise, which can cause significant thermal stress, and the total losses of the motor can be given by the following Equation [26]: cu al Fe + P P P P P P where, Pcu is winding loss, Pal is bar loss, PFe is core loss, P is mechanical loss, P is stray loss.Based on the results of electromagnetic field calculation, the losses can be obtained [27,28].The rotor core loss is relatively small and hence it can be ignored in Equation (8).

Results and Discussion
According to the above analysis, the coupled fluid-thermal field of the motor is calculated and the simulation results are analysed in this section.The absolute steady pressure based solver is adopted in this paper, and the solving time is about six hours.

Fluid Field Results
As can be seen from the Figure 8, the 3D fluid distribution on the shell surface, from the figure, the fluid continuous distribution can be seen.The air velocity on the shell surface appears to show a significant variation along the axis, and the max value of air velocity located in the fan side, it is about 10 m/s.The average air velocity on the shell surface is 4~5 m/s.Besides, the fluid status inside and outside of the junction box is obviously different under rated load, all of these explain why the motor temperature distribution along the axial and radial direction is not symmetrical.The air velocity along the path l1 was measured during the experiment, and the 19 testing points was selected (see Figure 9), and the test data are compared with the calculated results in the Figure 10.

Losses
As well known the losses, in an induction motor bring about the temperature rise, which can cause significant thermal stress, and the total losses of the motor can be given by the following Equation [26]: where, P cu is winding loss, P al is bar loss, P Fe is core loss, P Ω is mechanical loss, P ∆ is stray loss.
Based on the results of electromagnetic field calculation, the losses can be obtained [27,28].The rotor core loss is relatively small and hence it can be ignored in Equation (8).

Results and Discussion
According to the above analysis, the coupled fluid-thermal field of the motor is calculated and the simulation results are analysed in this section.The absolute steady pressure based solver is adopted in this paper, and the solving time is about six hours.

Fluid Field Results
As can be seen from the Figure 8, the 3D fluid distribution on the shell surface, from the figure, the fluid continuous distribution can be seen.The air velocity on the shell surface appears to show a significant variation along the axis, and the max value of air velocity located in the fan side, it is about 10 m/s.The average air velocity on the shell surface is 4~5 m/s.Besides, the fluid status inside and outside of the junction box is obviously different under rated load, all of these explain why the motor temperature distribution along the axial and radial direction is not symmetrical.The air velocity along the path l1 was measured during the experiment, and the 19 testing points was selected (see Figure 9), and the test data are compared with the calculated results in the Figure 10.   Figure 11 shows 2D streamlines of outer fluid field of the motor, from the figure, the air gets into fan cover and it will be affected by the rotating of the fan, then the air flows out along the heat radiation fins at high velocity, which leads to the forced heat exchange on shell surface of the motor.Figure 11 shows 2D streamlines of outer fluid field of the motor, from the figure, the air gets into fan cover and it will be affected by the rotating of the fan, then the air flows out along the heat radiation fins at high velocity, which leads to the forced heat exchange on shell surface of the motor.Figure 11 shows 2D streamlines of outer fluid field of the motor, from the figure, the air gets into fan cover and it will be affected by the rotating of the fan, then the air flows out along the heat radiation fins at high velocity, which leads to the forced heat exchange on shell surface of the motor.When fluid come to the load side of the motor, the velocity reduces rapidly due to the kinetic energy decrease, caused by open environment.Meanwhile, the obvious vortex phenomenon appears at region V1 and V2, which will lead to the energy loss of fluid medium and is one of the main reasons for vibration noise.Except for outer fluid field, the inner fluid field also is key to the temperature rise calculation.The TEFC induction motor is adopted in this paper, so the inner and outer fluid field are relatively separate, and the fluid flow inside the motor is completely caused by the rotor rotation.Figure 12 gives the 3D fluid distribution in the air-gap.The result shows that air velocity adjacent to rotor surface can reach about 10 m/s, and it is close to tangential speed of the rotor.Moreover, since stator is a stationary, and it has a strong effect on fluid movement in the air-gap.Therefore, the air velocity changes greatly, that is, it shows an obvious decline trend along the radial direction.
Energies 2018, 11, x FOR PEER REVIEW 9 of 17 appears at region V1 and V2, which will lead to the energy loss of fluid medium and is one of the main reasons for vibration noise.Except for outer fluid field, the inner fluid field also is key to the temperature rise calculation.The TEFC induction motor is adopted in this paper, so the inner and outer fluid field are relatively separate, and the fluid flow inside the motor is completely caused by the rotor rotation.Figure 12 gives the 3D fluid distribution in the air-gap.The result shows that air velocity adjacent to rotor surface can reach about 10 m/s, and it is close to tangential speed of the rotor.Moreover, since stator is a stationary, and it has a strong effect on fluid movement in the air-gap.Therefore, the air velocity changes greatly, that is, it shows an obvious decline trend along the radial direction.As shown in Figure 13, l2 is an axial path on the rotor surface, including the rotor core and the end-ring.The variation of air velocity along the path l2 is shown in Figure 14.From it, it is clear that the fluid field is generally symmetrical in the axial direction.The air speed on the rotor core surface has almost no change, and a substantial fluctuation can be found in the end-ring segment, which is attributed to violent turbulence in the confluence of the air gap and the end cavity.The accuracy of thermal analysis depends on the model.As a consequence, it also depends on the obtained coefficient of heat transfer related to internal and external heat exchange of the motor, including natural and forced convection.In this paper, a numerical algorithm of heat transfer coefficient is proposed, according to Newton's law of cooling and Fourier's law of heat conduction.As shown in Figure 13, l2 is an axial path on the rotor surface, including the rotor core and the end-ring.The variation of air velocity along the path l2 is shown in Figure 14.From it, it is clear that the fluid field is generally symmetrical in the axial direction.The air speed on the rotor core surface has almost no change, and a substantial fluctuation can be found in the end-ring segment, which is attributed to violent turbulence in the confluence of the air gap and the end cavity.
Energies 2018, 11, x FOR PEER REVIEW 9 of 17 appears at region V1 and V2, which will lead to the energy loss of fluid medium and is one of the main reasons for vibration noise.Except for outer fluid field, the inner fluid field also is key to the temperature rise calculation.The TEFC induction motor is adopted in this paper, so the inner and outer fluid field are relatively separate, and the fluid flow inside the motor is completely caused by the rotor rotation.Figure 12 gives the 3D fluid distribution in the air-gap.The result shows that air velocity adjacent to rotor surface can reach about 10 m/s, and it is close to tangential speed of the rotor.Moreover, since stator is a stationary, and it has a strong effect on fluid movement in the air-gap.Therefore, the air velocity changes greatly, that is, it shows an obvious decline trend along the radial direction.As shown in Figure 13, l2 is an axial path on the rotor surface, including the rotor core and the end-ring.The variation of air velocity along the path l2 is shown in Figure 14.From it, it is clear that the fluid field is generally symmetrical in the axial direction.The air speed on the rotor core surface has almost no change, and a substantial fluctuation can be found in the end-ring segment, which is attributed to violent turbulence in the confluence of the air gap and the end cavity.The accuracy of thermal analysis depends on the model.As a consequence, it also depends on the obtained coefficient of heat transfer related to internal and external heat exchange of the motor, including natural and forced convection.In this paper, a numerical algorithm of heat transfer coefficient is proposed, according to Newton's law of cooling and Fourier's law of heat conduction.appears at region V1 and V2, which will lead to the energy loss of fluid medium and is one of the main reasons for vibration noise.Except for outer fluid field, the inner fluid field also is key to the temperature rise calculation.The TEFC induction motor is adopted in this paper, so the inner and outer fluid field are relatively separate, and the fluid flow inside the motor is completely caused by the rotor rotation.Figure 12 gives the 3D fluid distribution in the air-gap.The result shows that air velocity adjacent to rotor surface can reach about 10 m/s, and it is close to tangential speed of the rotor.Moreover, since stator is a stationary, and it has a strong effect on fluid movement in the air-gap.Therefore, the air velocity changes greatly, that is, it shows an obvious decline trend along the radial direction.As shown in Figure 13, l2 is an axial path on the rotor surface, including the rotor core and the end-ring.The variation of air velocity along the path l2 is shown in Figure 14.From it, it is clear that the fluid field is generally symmetrical in the axial direction.The air speed on the rotor core surface has almost no change, and a substantial fluctuation can be found in the end-ring segment, which is attributed to violent turbulence in the confluence of the air gap and the end cavity.The accuracy of thermal analysis depends on the model.As a consequence, it also depends on the obtained coefficient of heat transfer related to internal and external heat exchange of the motor, including natural and forced convection.In this paper, a numerical algorithm of heat transfer coefficient is proposed, according to Newton's law of cooling and Fourier's law of heat conduction.The accuracy of thermal analysis depends on the model.As a consequence, it also depends on the obtained coefficient of heat transfer related to internal and external heat exchange of the motor, including natural and forced convection.In this paper, a numerical algorithm of heat transfer coefficient is proposed, according to Newton's law of cooling and Fourier's law of heat conduction.
After getting the initial fluid field, the initial heat transfer coefficient can be obtain, from which the heat field of the motor can be determined by the Equations ( 5)- (7).On account of the temperature of all nodes, a new heat transfer coefficient can be determined according to Equation ( 9), and the temperature of all nodes will be calculated again based on it: where i is the number of iterations, α i is coefficient of heat transfer for the ith iteration, λ is the thermal conductivity of cooling medium, T w is the surface temperature of solid located on the fluid-solid interface, T f is the average temperature of fluid, ∆X is the distance of two elements center, and it depends on the mesh generation strategies, T is temperature at position with distance of ∆X from the body surface.Equations ( 10) and ( 11) are used to identify whether the heat transfer coefficient and temperature are acceptable.The iteration is stopped if the heat transfer coefficients and temperatures meet the two convergence conditions of ( 10) and ( 11), otherwise the iteration continues, and the chart of the iterative method is shown in Figure 15.
Energies 2018, 11, x FOR PEER REVIEW 10 of 17 After getting the initial fluid field, the initial heat transfer coefficient can be obtain, from which the heat field of the motor can be determined by the Equations ( 5)- (7).On account of the temperature of all nodes, a new heat transfer coefficient can be determined according to Equation (9), and the temperature of all nodes will be calculated again based on it: where i is the number of iterations, αi is coefficient of heat transfer for the ith iteration, λ is the thermal conductivity of cooling medium, Tw is the surface temperature of solid located on the fluid-solid interface, Tf is the average temperature of fluid, ΔX is the distance of two elements center, and it depends on the mesh generation strategies, T is temperature at position with distance of ΔX from the body surface.Equations ( 10) and ( 11) are used to identify whether the heat transfer coefficient and temperature are acceptable.The iteration is stopped if the heat transfer coefficients and temperatures meet the two convergence conditions of ( 10) and ( 11), otherwise the iteration continues, and the chart of the iterative method is shown in Figure 15.( ) (α i+1 − α i ) where T k is the temperature obtained at the k times, ε 1 and ε 2 are residual error.To get a more precision heat transfer coefficient, the environmental temperature was chosen as the initial temperature of fluid field.The heat transfer coefficient distribution of the motor is displayed in Figure 16a,b.
The calculated results of fluid field can provide accurate boundary conditions for the thermal field calculation.From the results, it is revealed that the fluid field distribution has little change before and after the broken bar fault, so only the results of fluid field about healthy motor were reported.
Energies 2018, 11, x FOR PEER REVIEW 11 of 17 where T k is the temperature obtained at the k times, ε1 and ε2 are residual error.To get a more precision heat transfer coefficient, the environmental temperature was chosen as the initial temperature of fluid field.The heat transfer coefficient distribution of the motor is displayed in Figure 16a,b.The calculated results of fluid field can provide accurate boundary conditions for the thermal field calculation.From the results, it is revealed that the fluid field distribution has little change before and after the broken bar fault, so only the results of fluid field about healthy motor were reported.

Thermal Field Results
In this section, 3D steady thermal field of the motor operating at the full load is studied.Figure 17a-c shows the temperature distribution of the motor with healthy cage, one and two broken bars respectively.From the figures, the thermal field of the motor can be divided into two regions, the rotor part and the stator part, and the rotor temperature is much higher than the stator temperature.Because of the existence of the fan and the junction box, the temperature distribution of the motor shows obvious asymmetry along the axial and circumferential direction.Due to the good cooling conditions, the temperature of motor frame is the lowest, relative to the other parts.

Thermal Field Results
In this section, 3D steady thermal field of the motor operating at the full load is studied.Figure 17a-c shows the temperature distribution of the motor with healthy cage, one and two broken bars respectively.From the figures, the thermal field of the motor can be divided into two regions, the rotor part and the stator part, and the rotor temperature is much higher than the stator temperature.Because of the existence of the fan and the junction box, the temperature distribution of the motor shows obvious asymmetry along the axial and circumferential direction.Due to the good cooling conditions, the temperature of motor frame is the lowest, relative to the other parts.
Energies 2018, 11, x FOR PEER REVIEW 11 of 17 where T k is the temperature obtained at the k times, ε1 and ε2 are residual error.To get a more precision heat transfer coefficient, the environmental temperature was chosen as the initial temperature of fluid field.The heat transfer coefficient distribution of the motor is displayed in Figure 16a,b.The calculated results of fluid field can provide accurate boundary conditions for the thermal field calculation.From the results, it is revealed that the fluid field distribution has little change before and after the broken bar fault, so only the results of fluid field about healthy motor were reported.

Thermal Field Results
In this section, 3D steady thermal field of the motor operating at the full load is studied.Figure 17a-c shows the temperature distribution of the motor with healthy cage, one and two broken bars respectively.From the figures, the thermal field of the motor can be divided into two regions, the rotor part and the stator part, and the rotor temperature is much higher than the stator temperature.Because of the existence of the fan and the junction box, the temperature distribution of the motor shows obvious asymmetry along the axial and circumferential direction.Due to the good cooling conditions, the temperature of motor frame is the lowest, relative to the other parts.The rotary magnetic field is created by the stator current, which brings a power loss at the same time, and stator winding temperature rise caused by the loss is always concerned.Figure 18 gives the temperature of the stator winding in different cases.The highest temperature of the stator winding with healthy cage is about 76 • C and it increases obviously in case of faulty condition which will result in insulation aging.Besides, because of the asymmetric cooling condition on the shell surface, the temperature of the stator winding close to the junction box is higher than other windings.The rotary magnetic field is created by the stator current, which brings a power loss at the same time, and stator winding temperature rise caused by the loss is always concerned.Figure 18 gives the temperature of the stator winding in different cases.The highest temperature of the stator winding with healthy cage is about 76 °C and it increases obviously in case of faulty condition which will result in insulation aging.Besides, because of the asymmetric cooling condition on the shell surface, the temperature of the stator winding close to the junction box is higher than other windings.The temperature distributions of the rotor with healthy bars and broken bars are given in Figure 19.Compared with other parts of the motor, the temperature is the highest.Due to the good thermal conductivities of the bars and rotor core, there is minimal temperature changes of the entire rotor can be observed.The temperature distributions of the rotor with healthy bars and broken bars are given in Figure 19.Compared with other parts of the motor, the rotor temperature is the highest.Due to the good thermal conductivities of the bars and rotor core, there is minimal temperature changes of the entire rotor can be observed.The rotary magnetic field is created by the stator current, which brings a power loss at the same time, and stator winding temperature rise caused by the loss is always concerned.Figure 18 gives the temperature of the stator winding in different cases.The highest temperature of the stator winding with healthy cage is about 76 °C and it increases obviously in case of faulty condition which will result in insulation aging.Besides, because of the asymmetric cooling condition on the shell surface, the temperature of the stator winding close to the junction box is higher than other windings.The temperature distributions of the rotor with healthy bars and broken bars are given in Figure 19.Compared with other parts of the motor, the rotor temperature is the highest.Due to the good thermal conductivities of the bars and rotor core, there is minimal temperature changes of the entire rotor can be observed.It can be seen from the whole rotor solving region that the lowest temperatures are located at broken bars, as shown in Figure 19b,c.Despite there is no currents pass through the broken bars, as well as no aluminium losses are produced.Based on the magnetic field calculation results [27][28][29], it is revealed that the current and the aluminium losses in the rotor bar adjacent to the faulty bars increase obviously.The current distribution of rotor bars for different cases is shown in Figure 20a.But in fact, there are little temperature differences between the healthy bars and broken bars, even the whole rotor, due to the good thermal conductivity properties of the rotor materials.It can be seen from the whole rotor solving region that the lowest temperatures are located at broken bars, as shown in Figure 19b,c.Despite there is no currents pass through the broken bars, as well as no aluminium losses are produced.Based on the magnetic field calculation results [27][28][29], it is revealed that the current and the aluminium losses in the rotor bar adjacent to the faulty bars increase obviously.The current distribution of rotor bars for different cases is shown in Figure 20a.But in fact, there are little temperature differences between the healthy bars and broken bars, even the whole rotor, due to the good thermal conductivity properties of the rotor materials.The temperature of some representative measured points is listed in Table 4, and the test values and calculated values are compared.Because it took a long time to do the temperature experiments, the change of ambient temperature must to be considered.In this paper, the ambient temperature at above three cases are 20.5 °C, 16.5 °C, 17 °C respectively.From the data in Table 4, it is cleared that the temperature rise is closely related to the broken bar number, that is, it increases rapidly with increase of broken bar number.In addition, the calculated values of temperature rise are similar to the test values at measurement point, and the accuracy of the simulation results is proved.The temperature of some representative measured points is listed in Table 4, and the test values and calculated values are compared.Because it took a long time to do the temperature experiments, the change of ambient temperature must to be considered.In this paper, the ambient temperature at above three cases are 20.5 • C, 16.5 • C, 17 • C respectively.From the data in Table 4, it is cleared that the temperature rise is closely related to the broken bar number, that is, it increases rapidly with increase of broken bar number.In addition, the calculated values of temperature rise are similar to the test values at measurement point, and the accuracy of the simulation results is proved.Figure 21 shows the temperature rise distributions of the motor along the radial direction for the three cases described above.It is indicates that the change of the temperature rise is visible.Except for the temperature increase due to the fault, there is subtle difference in radial temperature distribution tendencies of the motor for the three different cases.There is a large temperature difference between the rotor and stator owing to the poor thermal conductivity of air, and the same situation also exists between the stator and motor shell because of the contact thermal resistance.
Figure 21 shows the temperature rise distributions of the motor along the radial direction for the three cases described above.It is indicates that the change of the temperature rise is visible.Except for the temperature increase due to the fault, there is subtle difference in radial temperature distribution tendencies of the motor for the three different cases.There is a large temperature difference between the rotor and stator owing to the poor thermal conductivity of air, and the same situation also exists between the stator and motor shell because of the contact thermal resistance.

Conclusions
In this paper, the research developed a CFD model for coupled fluid-thermal calculation of the medium small induction motor with healthy and broken bar fault rotors.The fluid field and the thermal field have been analyzed around the squirrel cage with or without the broken bar fault.Compared with the healthy case, the growth of these quantities is not obvious in the one broken-bar case but they increase greatly in two adjacent broken bars, and the thermal field distributions are asymmertrical in the condition of fault cases.The calculated results of temperature regarding the  Figure 21 shows the temperature rise distributions of the motor along the radial direction for the three cases described above.It is indicates that the change of the temperature rise is visible.Except for the temperature increase due to the fault, there is subtle difference in radial temperature distribution tendencies of the motor for the three different cases.There is a large temperature difference between the rotor and stator owing to the poor thermal conductivity of air, and the same situation also exists between the stator and motor shell because of the contact thermal resistance.

Conclusions
In this paper, the research developed a CFD model for coupled fluid-thermal calculation of the medium small induction motor with healthy and broken bar fault rotors.The fluid field and the thermal field have been analyzed around the squirrel cage with or without the broken bar fault.Compared with the healthy case, the growth of these quantities is not obvious in the one broken-bar case but they increase greatly in two adjacent broken bars, and the thermal field distributions are asymmertrical in the condition of fault cases.The calculated results of temperature regarding the

Conclusions
In this paper, the research developed a CFD model for coupled fluid-thermal calculation of the medium small induction motor with healthy and broken bar fault rotors.The fluid field and the thermal field have been analyzed around the squirrel cage with or without the broken bar fault.Compared with the healthy case, the growth of these quantities is not obvious in the one broken-bar case but they increase greatly in two adjacent broken bars, and the thermal field distributions are asymmertrical in the condition of fault cases.The calculated results of temperature regarding the fluid field are compared with the test results, shows the accuracy and rationality of the model.Furthermore,

Figure 1 .
Figure 1.The laboratory test setup and the broken bar rotors: (a) The laboratory test setup; (b) the broken bar rotors.

Figure 1 .
Figure 1.The laboratory test setup and the broken bar rotors: (a) The laboratory test setup; (b) the broken bar rotors.

Figure 1 .
Figure 1.The laboratory test setup and the broken bar rotors: (a) The laboratory test setup; (b) the broken bar rotors.

Figure 2 .
Figure 2. 3D model of the motor.Figure 2. 3D model of the motor.

Figure 2 .
Figure 2. 3D model of the motor.Figure 2. 3D model of the motor.

Figure 3 .
Figure 3.The sketch map of the equivalent stator slots.

Figure 4 .
Figure 4.The structure in the modeling of the equivalent stator slots

Figure 3 .
Figure 3.The sketch map of the equivalent stator slots.

Figure 3 .
Figure 3.The sketch map of the equivalent stator slots.

Figure 4 .
Figure 4.The structure in the modeling of the equivalent stator slots

Figure 4 .
Figure 4.The structure in the modeling of the equivalent stator slots.
generation ratio, and σ k and σ e are turbulence Prandtl number (Pr number) of the k equation and the ε equation, respectively.

Figure 5 .
Figure 5.The model of external fluid domain of the motor.

Figure 6 .
Figure 6.Rotating fluid domain of fan.

Figure 7
Figure 7 is the section diagram of air-gap fluid field, and it is divided into two regions, including the rotary and stationary fluid region.The rotary fluid region is close to the rotor, and the stationary fluid region is close to the stator.On the basis of the reference coordinate system, the

Figure 5 .
Figure 5.The model of external fluid domain of the motor.

Figure 5 .
Figure 5.The model of external fluid domain of the motor.

Figure 6 .
Figure 6.Rotating fluid domain of fan.

Figure 7
Figure 7 is the section diagram of air-gap fluid field, and it is divided into two regions, including the rotary and stationary fluid region.The rotary fluid region is close to the rotor, and the stationary fluid region is close to the stator.On the basis of the reference coordinate system, the

Figure 6 .
Figure 6.Rotating fluid domain of fan.

Figure 7 .
Figure 7.The section diagram of air-gap fluid field.

Figure 7 .Table 3 .
Figure 7.The section diagram of air-gap fluid field.

Figure 8 .
Figure 8. 3D fluid field distribution on the shell surface of the motor.

Figure 9 .
Figure 9. Air velocity measured positions on the motor frame.

Figure 10 .
Figure 10.Comparison of air velocity in simulation and experiment along the path 11.

Figure 11 .
Figure 11.2D streamlines of outer fluid field of the motor.When fluid come to the load side of the motor, the velocity reduces rapidly due to the kinetic energy decrease, caused by open environment.Meanwhile, the obvious vortex phenomenon

Figure 8 . 17 Figure 8 .
Figure 8. 3D fluid field distribution on the shell surface of the motor.

l 1 1 19Figure 9 .
Figure 9. Air velocity measured positions on the motor frame.

Figure 10 .
Figure 10.Comparison of air velocity in simulation and experiment along the path 11.

Figure 11 .
Figure 11.2D streamlines of outer fluid field of the motor.When fluid come to the load side of the motor, the velocity reduces rapidly due to the kinetic energy decrease, caused by open environment.Meanwhile, the obvious vortex phenomenon

Figure 9 . 17 Figure 8 .
Figure 9. Air velocity measured positions on the motor frame.

l 1 1 19Figure 9 .
Figure 9. Air velocity measured positions on the motor frame.

Figure 10 .
Figure 10.Comparison of air velocity in simulation and experiment along the path 11.

Figure 11 .
Figure 11.2D streamlines of outer fluid field of the motor.When fluid come to the load side of the motor, the velocity reduces rapidly due to the kinetic energy decrease, caused by open environment.Meanwhile, the obvious vortex phenomenon

Figure 10 .
Figure 10.Comparison of air velocity in simulation and experiment along the path 11.

Figure 11 17 Figure 8 .
Figure 11  shows 2D streamlines of outer fluid field of the motor, from the figure, the air gets into fan cover and it will be affected by the rotating of the fan, then the air flows out along the heat radiation fins at high velocity, which leads to the forced heat exchange on shell surface of the motor.

l 1 1 19Figure 9 .
Figure 9. Air velocity measured positions on the motor frame.

Figure 10 .
Figure 10.Comparison of air velocity in simulation and experiment along the path 11.

Figure 11 shows
Figure11shows 2D streamlines of outer fluid field of the motor, from the figure, the air gets into fan cover and it will be affected by the rotating of the fan, then the air flows out along the heat radiation fins at high velocity, which leads to the forced heat exchange on shell surface of the motor.

Figure 11 .
Figure 11.2D streamlines of outer fluid field of the motor.

Figure 11 .
Figure 11.2D streamlines of outer fluid field of the motor.

Figure 12 .
Figure 12. 3D fluid field distribution in the air-gap.

Figure 13 .
Figure 13.3D model of the rotor.

Figure 14 .
Figure 14.Fluid velocity distribution on rotor surface along axial.

Figure 12 .
Figure 12. 3D fluid field distribution in the air-gap.

Figure 12 .
Figure 12. 3D fluid field distribution in the air-gap.

Figure 13 .
Figure 13.3D model of the rotor.

Figure 14 .
Figure 14.Fluid velocity distribution on rotor surface along axial.

Figure 13 .
Figure 13.3D model of the rotor.

Figure 12 .
Figure 12. 3D fluid field distribution in the air-gap.

Figure 13 .
Figure 13.3D model of the rotor.

Figure 14 .
Figure 14.Fluid velocity distribution on rotor surface along axial.

Figure 14 .
Figure 14.Fluid velocity distribution on rotor surface along axial.

Figure 15 .
Figure 15.The chart of the iterative method.

Figure 15 .
Figure 15.The chart of the iterative method.

Figure 16 .
Figure 16.Heat transfer coefficient distribution of the motor: (a) Heat transfer surfaces distribution of the stator and frame; (b) Heat transfer coefficient distribution of the rotor.

Figure 16 .
Figure 16.Heat transfer coefficient distribution of the motor: (a) Heat transfer surfaces distribution of the stator and frame; (b) Heat transfer coefficient distribution of the rotor.

Figure 16 .
Figure 16.Heat transfer coefficient distribution of the motor: (a) Heat transfer surfaces distribution of the stator and frame; (b) Heat transfer coefficient distribution of the rotor.

Figure 20 .
Figure 20.The current distribution: (a) Healthy motor cage and one-broken-bar fault; (b) Healthy motor cage and continuous two-broken-bars.

Figure 20 .
Figure 20.The current distribution: (a) Healthy motor cage and one-broken-bar fault; (b) Healthy motor cage and continuous two-broken-bars.

Figure 21 .
Figure 21.Temperature rise distributions of the motor along the radial direction for the three cases.

Figure 22 Figure 22 .
Figure 22 displays the heat flux distribution of the rotor ring and bar, from which, the heat flux distribution of the rotor appears to be uneven, and heat flux value of the ring is slightly larger than that of rotor bars.There is fairly high heat flux in the junction of the bars and the rings, where the broken bar fault frequently occurs.Meanwhile, when broken bar fault happens, the heat flux value of the joint adjacent to fault bars goes up greatly, which will probably result in aggravation of the broken bar fault.

Figure 21 .
Figure 21.Temperature rise distributions of the motor along the radial direction for the three cases.

Figure 22
Figure 22 displays the heat flux distribution of the rotor ring and bar, from which, the heat flux distribution of the rotor appears to be uneven, and heat flux value of the ring is slightly larger than that of rotor bars.There is fairly high heat flux in the junction of the bars and the rings, where the broken bar fault frequently occurs.Meanwhile, when broken bar fault happens, the heat flux value of the joint adjacent to fault bars goes up greatly, which will probably result in aggravation of the broken bar fault.

Figure 21 .
Figure 21.Temperature rise distributions of the motor along the radial direction for the three cases.

Figure 22 Figure 22 .
Figure 22 displays the heat flux distribution of the rotor ring and bar, from which, the heat flux distribution of the rotor appears to be uneven, and heat flux value of the ring is slightly larger than that of rotor bars.There is fairly high heat flux in the junction of the bars and the rings, where the broken bar fault frequently occurs.Meanwhile, when broken bar fault happens, the heat flux value of the joint adjacent to fault bars goes up greatly, which will probably result in aggravation of the broken bar fault.

Table 1 .
Basic parameters of the prototype motor.

Table 1 .
Basic parameters of the prototype motor.

Table 1 .
Basic parameters of the prototype motor.

Table 2 .
Properties of material.

Table 3 .
Boundary conditions of computational fluid dynamic (CFD) simulation.

Table 4 .
Comparison of the temperature value of the measured points (°C).

Table 4 .
Comparison of the temperature value of the measured points ( • C).
eq equivalent heat conductivity of the insulation layer in the stator slots, W/m•k