An Approximate Estimation Approach of Fault Size for Spalled Ball Bearing in Induction Motor by Tracking Multiple Vibration Frequencies in Current

Fault size estimation is of great importance to bearing performance degradation assessment and life prediction. Until now, fault size estimation has generally been based on acoustic emission signals or vibration signals; an approach based on current signals has not yet been mentioned. In the present research, an approximate estimation approach based on current is introduced. The proposed approach is easy to implement for existing inverter-driven induction motors without complicated calculations and additional sensors, immune to external disturbances, and suitable for harsh conditions. Firstly, a feature transmission route from spall, to Hertzian forces, and then to friction torque is simulated based on a spall model and dynamic model of the bearing. Based on simulated results, the relation between spall size and the multiple characteristic vibration frequencies in friction torque is revealed. Secondly, the multiple characteristic vibration frequencies modulated in the current is investigated. Analysis results show that those frequencies modulated in the current are independent of each other, without spectrum overlap. Thirdly, to address the issue of which fault features modulated in the current are very weak, a fault-feature-highlighting approach based on reduced voltage frequency ratio is proposed. Finally, experimental tests were conducted. The obtained results validate that the proposed approach is feasible and effective for spall size estimation.


Introduction
Induction motors play an important role in many industrial and on-board applications because of their low cost, simple construction, and high reliability. It is well known that faulty bearings contribute to most of the failures in rotating machinery [1,2]. Depending on the type and size of the machine, bearing failure distributions vary from about 40% to about 90% from large to small machines [3]. In fact, bearings, even when properly designed, are sensitive components and failure is often due to inadequate operating conditions or failures in maintenance, which are more serious on marine ships than on land, such as excessive loading, shaft misalignment, wrong mounting, improper lubrication, etc. [4].
Fatigue in rolling element bearings, resulting in spalling of the races and rolling elements, is the most common cause of bearing failures. There are three types of fatigue in bearings: surface distress, fatigue pitting and fatigue spalling [5]. Surface distress appears as a smooth surface resulting from plastic deformation in the asperity region (typically less than 10 µm). Pitting appears as shallow craters in contact surfaces with a depth of, at most, the thickness of the work-hardened layer (approximately 10 µm). Spalling leaves deeper cavities at contact surfaces with a depth of 20-100 µm [6].
Bearing faults are generally slowly progressive [7]. The appearance of the first spall does not cause immediate breakdown [8], and also does not mean the end of the bearing's useful life [9]. A premature removal of the bearing from service may be very expensive [6]. The size of the spalling area has a significant influence on the operation performance and the remaining useful life of rolling element bearings. Therefore, spall size estimation is of great importance to bearing performance degradation assessment and life prediction [9]. If the spall size can be promptly tracked, it is greatly helpful for ensuring correct maintenance decisions and preventing unexpected downtime [10].
For the moment, the spall size estimation approach is generally based on acoustic emission signal [2,9,11], or vibration signal [6,10,[12][13][14]. Kang et al. [2] introduced a 2-D visualization tool that represents the percentage of the Gaussian-mixture-model-based residual component-to-defect component ratios via time-varying and multi-resolution envelope analysis, then distinguished the different crack sizes based on the k-NN classifier, but the specific fault size was not estimated. Ming et al. [9] used an averaged dual-impulse interval determining method to evaluate the spall size by calculating the autocorrelation function of the squared envelope. Al-Ghamd et al. [11] used the burst duration in acoustic emission signals to estimate the spall size. Sawalhi et al. [6] used cepstrum analysis to give an average estimate of the spacing between the entry and impact events. Zhao et al. [10] used the approximate entropy method and empirical mode decomposition to extract the times of the entry and exit events. Cui et al. [12] used the vertical-horizontal synchronized root mean square to distinguish different fault sizes. Ahmadi et al. [13,14] estimated the defect size by measuring the distance between the entry and impact events. Most of the above approaches implement fault size estimation by determining the interval between two events caused by the same spall.
A common drawback of the above approaches is that they are not suitable for applications in some special and harsh environments which have external vibration or noise interference. In these cases, current signals are preferable as they are immune to external environmental interference. However, the two events caused by the spall are modulated in the alternating current, and are scarcely possible to distinguished in the time domain using the above approaches. Therefore, spall size estimation based on current has not been mentioned until now.
A novel idea which tracks the two events in the frequency domain is proposed in the present work. To address this issue, the feature transmission route from spall, to Hertzian forces, to friction torque, and then to current, is investigated detailedly. Then, the relation between the spall size and the sidebands modulated in current is revealed. However, identifying those sidebands modulated in current is consistently very tricky.
The complex signal transmission route from bearing defect to stator current poses higher impedance and hence results in lower signal-to-noise ratio (SNR) [15]. The amplitude of the torque ripple is very small for real defects; the coefficient between the torque ripple and the current modulation is a damping coefficient [3]. The modulation index is largely less than one [16]. Owing to the poor SNR, the detection and diagnosis of bearings faults by using current signal is still challenging [7].
Furthermore, the current of a typical induction motor involves dominant components such as the supply fundamental and its multiple harmonics, the eccentricity, slot, and saturation harmonics, etc. [17]. For the defect signatures, the dominant components in current are not related to bearing defects, which are usually considered as noise. The suppression of strong noise is then another problem to be solved. Therefore, highlighting fault features (enhancing SNR) or noise suppression are always a research focus, whether in current or in other signals. Noise suppression in current-based bearing fault diagnosis is focused on elimination of the fundamental supply frequency and its harmonics. To address this problem, some effective approaches are proposed based on notch filters [18], the noise cancellation method using time shifting [17], and the Teager-Kaiser energy operator [16], etc. However, as the power supply frequency is the carrier of defect signatures, excessive noise suppression is not always the most reasonable option. Instead, a simple and effective approach of noise suppression is the squared envelope of signal based on Hilbert transform. The squared envelope of signal retains all the necessary Sensors 2020, 20, 1631 3 of 23 diagnostic information [16], and partly suppresses the power supply frequency and its harmonics. Till now, the squared envelope of signal has been applied successfully in preprocessing of current and other signals [3,17,[19][20][21].
Focusing on the poor SNR, some effective approaches to highlighting fault features in current have been proposed based on space vector [3], Park's vector [22], etc. Unlike the above conventional strategies, this paper attempts to highlight fault features based on motor operation under reduced voltage frequency ratio.
Finally, a simple and feasible spall size estimation approach by tracking multiple vibration frequencies in current is introduced. The proposed approach only needs to calculate the squared envelope spectrum, an algorithm which is very simple, mature and widely used in various fields. Therefore, the proposed approach is easy to implement for existing inverter-driven induction motors without complicated calculation and additional sensors, immune to external disturbances, and suitable for harsh conditions. Two major contributions have been made in this work. First, the relation between the spall size and the sidebands modulated in the current is revealed based on the simulation results for a spall model and a dynamic model of bearings; then, a novel idea which estimates the spall size by tracking the current sidebands is introduced. Second, this work introduces a fault-feature-highlighting approach based on reduced voltage frequency ratio, which is more effective than the traditional approaches.
This work is organized as follows. In Section 2, the fluctuation of friction torque caused by spall, the relation between spall size and multiple characteristic vibration frequencies, and the multiple characteristic vibration frequencies modulated in current, are investigated. In Section 3, the fluctuation of friction torque is verified based on a customized measuring instrument, and the spall size estimation approach by tracking multiple vibration frequencies is verified based on experimental platform. Finally, the conclusions and remarks are given in Section 4.

Elastic Deformation
For ball bearings, widely employed in induction motors, elastic deformation between raceways and balls produces a non-linear phenomenon between force and deformation, which is obtained by the Hertzian theory [23]. The non-linear relation of load-deformation is given by where F is Hertzian contact force, K is the load-deflection factor or constant for Hertzian contact elastic deformation, δ is the radial deflection or contact deformation. The load-deflection factor K depends on the contact geometry. Total deflection between two raceways is the sum of the approaches between the rolling elements and each raceway. K is given by [24].
where K i is inner-raceway-to-ball contact stiffness, and K o is outer-raceway-to-ball contact stiffness. The detailed formulas and parameters for K i and K o are introduced in the literature [23][24][25]. The value of K for an SKF 6206 bearing is 4.9582 × 107 N/m, according with to the value range proposed in [25][26][27].

Spall Model of Bearings
We considered the spall in the outer raceway as a rectangular notch of depth h d with a width of w d and a length of L d shown in Figure 1. Take into account that the spall depth is generally larger than the radial deflection of outer raceways, and the spall width is also larger than the diameter of Sensors 2020, 20, 1631 4 of 23 the contact deformation area between raceways and balls. Therefore, the spall size estimation in the present research focused on different spall lengths, L d .
Sensors 2020, 20, 1631 4 of 23 the contact deformation area between raceways and balls. Therefore, the spall size estimation in the present research focused on different spall lengths, Ld.
(a) (b) The rolling element-raceway contact can be considered as a spring mass system, in which the outer race is fixed in a rigid support and the inner race is fixed rigidly with the motor shaft. In Figure  1a, the 7th ball is shown as the symbol of spring mass, the other balls are only shown as a number for convenience.
Generally, ball bearings are designed with clearance Cr for flexible and unblocked operation [23]. On account of the effect of clearance Cr, the number of balls providing rotor support varies. If Cr = 0, there are five balls (4,3,2,1,9) providing rotor support as shown in Figure 1a. If Cr > 0, there are three balls (3, 2, 1) providing rotor support. For an SKF 6206 bearing, the clearance Cr is generally about 5-20 μm. Therefore, the balls in the area between 5 o'clock and 7 o'clock principally provide rotor support; it is clear that the effect on the motor is most severe for the spall located in this support area.
Because the load-deflection factor K is very large, the contact deformation δ is correspondingly very small, and generally micrometer-scale. The contact deformation area between the ball and the race is very small; the rise and fall of contact force will be obvious even for a very small spall. When the 2nd ball rolls into the spall, the contact force F2 of the 2nd ball disappears, the centre of the inner race connected with the motor shaft falls, the contact force F1 and F3 subsequently rises, and the vibration is also produced at this moment.
For further exploration of the effect of spall, a detailed rolling path is shown in Figure 2. Considering the balls as point masses is a common approach employed in numerous literatures [24,[28][29][30], in which the contact force of the ball entirely disappears as soon as the centre of the ball enters into the entry of the spall. When the ball passes through the spall from the entry to the exit of the spall, the common rolling path of the centre of the ball is considered to be in a sequence of "a-bc-d-e-f-g-h" as shown in Figure 2. In reality, considering the size of the ball [25,[31][32][33][34], the rolling path of the centre of the ball is in a sequence of "a-b-d-e-g-h". Till the "d" point, the contact force of the ball just entirely disappears; beginning at the "e" point, the contact force of the ball already arises. The rolling element-raceway contact can be considered as a spring mass system, in which the outer race is fixed in a rigid support and the inner race is fixed rigidly with the motor shaft. In Figure 1a, the 7th ball is shown as the symbol of spring mass, the other balls are only shown as a number for convenience.
Generally, ball bearings are designed with clearance C r for flexible and unblocked operation [23]. On account of the effect of clearance C r , the number of balls providing rotor support varies. If C r = 0, there are five balls (4,3,2,1,9) providing rotor support as shown in Figure 1a. If C r > 0, there are three balls (3, 2, 1) providing rotor support. For an SKF 6206 bearing, the clearance C r is generally about 5-20 µm. Therefore, the balls in the area between 5 o'clock and 7 o'clock principally provide rotor support; it is clear that the effect on the motor is most severe for the spall located in this support area.
Because the load-deflection factor K is very large, the contact deformation δ is correspondingly very small, and generally micrometer-scale. The contact deformation area between the ball and the race is very small; the rise and fall of contact force will be obvious even for a very small spall. When the 2nd ball rolls into the spall, the contact force F 2 of the 2nd ball disappears, the centre of the inner race connected with the motor shaft falls, the contact force F 1 and F 3 subsequently rises, and the vibration is also produced at this moment.
For further exploration of the effect of spall, a detailed rolling path is shown in Figure 2.
Sensors 2020, 20, 1631 4 of 23 the contact deformation area between raceways and balls. Therefore, the spall size estimation in the present research focused on different spall lengths, Ld.
(a) (b) The rolling element-raceway contact can be considered as a spring mass system, in which the outer race is fixed in a rigid support and the inner race is fixed rigidly with the motor shaft. In Figure  1a, the 7th ball is shown as the symbol of spring mass, the other balls are only shown as a number for convenience.
Generally, ball bearings are designed with clearance Cr for flexible and unblocked operation [23]. On account of the effect of clearance Cr, the number of balls providing rotor support varies. If Cr = 0, there are five balls (4,3,2,1,9) providing rotor support as shown in Figure 1a. If Cr > 0, there are three balls (3, 2, 1) providing rotor support. For an SKF 6206 bearing, the clearance Cr is generally about 5-20 μm. Therefore, the balls in the area between 5 o'clock and 7 o'clock principally provide rotor support; it is clear that the effect on the motor is most severe for the spall located in this support area.
Because the load-deflection factor K is very large, the contact deformation δ is correspondingly very small, and generally micrometer-scale. The contact deformation area between the ball and the race is very small; the rise and fall of contact force will be obvious even for a very small spall. When the 2nd ball rolls into the spall, the contact force F2 of the 2nd ball disappears, the centre of the inner race connected with the motor shaft falls, the contact force F1 and F3 subsequently rises, and the vibration is also produced at this moment.
For further exploration of the effect of spall, a detailed rolling path is shown in Figure 2. Considering the balls as point masses is a common approach employed in numerous literatures [24,[28][29][30], in which the contact force of the ball entirely disappears as soon as the centre of the ball enters into the entry of the spall. When the ball passes through the spall from the entry to the exit of the spall, the common rolling path of the centre of the ball is considered to be in a sequence of "a-bc-d-e-f-g-h" as shown in Figure 2. In reality, considering the size of the ball [25,[31][32][33][34], the rolling path of the centre of the ball is in a sequence of "a-b-d-e-g-h". Till the "d" point, the contact force of the ball just entirely disappears; beginning at the "e" point, the contact force of the ball already arises. Considering the balls as point masses is a common approach employed in numerous literatures [24,[28][29][30], in which the contact force of the ball entirely disappears as soon as the centre of the ball enters into the entry of the spall. When the ball passes through the spall from the entry to the exit of the spall, the common rolling path of the centre of the ball is considered to be in a sequence of "a-b-c-d-e-f -g-h" as shown in Figure 2. In reality, considering the size of the ball [25,[31][32][33][34], the rolling path of the centre of the ball is in a sequence of "a-b-d-e-g-h". Till the "d" point, the contact force of the ball just entirely disappears; beginning at the "e" point, the contact force of the ball already arises. In order to estimate accurately the effect of spall, the present research adopted the spall model introduced by Alireza Moazen Ahmadi [35], which reasonably considers the size of the ball. The schematic diagram of the ball passing the spall in the outer raceway is as shown in Figure 3. In order to estimate accurately the effect of spall, the present research adopted the spall model introduced by Alireza Moazen Ahmadi [35], which reasonably considers the size of the ball. The schematic diagram of the ball passing the spall in the outer raceway is as shown in Figure 3. The spall in the outer raceway can be modeled as where θ is the angular position of the ball, θen and θex are the angular positions of the spall entry and exit, and hd is the spall depth. The geometry function of the outer raceway is given by where rc is bearing pitch radius, rb is the radius of the ball, and Cr is the radial clearance of the bearing. The contact deformation δ is given by where Z is the distance of the ball from the centre of the outer raceway, β is the angle between the maximum deformation point on the ball and Z based on the centre of the ball, and ψ is the angle between the maximum deformation point on the ball and Z based on the centre of the outer raceway. The ψ is given by The contact deformation δd between the entry and the exit of the spall is finally determined by

Dynamic Model of Bearings
Because the ball bearing primarily supports the radial forces, the dynamic model based on a two-degrees-of-freedom system is enough for analysing the defect feature transmission route from spall, to Hertzian forces, and then to friction torque. The two-degrees-of-freedom system is also employed widely in the majority of literature up to now [24][25][28][29][30][31][32]. The dynamic model of bearings is shown in Figure 4. Φd is the angle from the X-axis to the centre of the spall, and ΔΦd is the angle corresponding to the spall length Ld. For convenience, ΔΦd is directly regarded as the spall length in the following part. The spall in the outer raceway can be modeled as where θ is the angular position of the ball, θ en and θ ex are the angular positions of the spall entry and exit, and h d is the spall depth. The geometry function of the outer raceway is given by where r c is bearing pitch radius, r b is the radius of the ball, and C r is the radial clearance of the bearing. The contact deformation δ is given by where Z is the distance of the ball from the centre of the outer raceway, β is the angle between the maximum deformation point on the ball and Z based on the centre of the ball, and ψ is the angle between the maximum deformation point on the ball and Z based on the centre of the outer raceway. The ψ is given by The contact deformation δ d between the entry and the exit of the spall is finally determined by

Dynamic Model of Bearings
Because the ball bearing primarily supports the radial forces, the dynamic model based on a two-degrees-of-freedom system is enough for analysing the defect feature transmission route from spall, to Hertzian forces, and then to friction torque. The two-degrees-of-freedom system is also employed widely in the majority of literature up to now [24,25,[28][29][30][31][32]. The dynamic model of bearings is shown in Figure 4. Φ d is the angle from the X-axis to the centre of the spall, and ∆Φ d is the angle If x and y are the deflections along the X-axis and the Y-axis, the radial deflection δi of the i th ball at any angle θi is given by The Z in Equation (5) (10) where N is the number of balls. Since the Hertzian contact force arises only when there is contact between the ball and the raceway, the respective contact force is set to zero when the contact deformation is equal or smaller than zero. This is indicated by subscript "+". θi is given by where ωc is the speed of the cage.
where ω is the shaft speed. Based on the structural dynamics theory, the equations of motion for a two-degrees-of-freedom system can be written as follows:  If x and y are the deflections along the X-axis and the Y-axis, the radial deflection δ i of the i th ball at any angle θ i is given by The Z in Equation (5) is given by Based on Equation (1), resolving the total Hertzian forces along the X-axis and Y-axis: where N is the number of balls. Since the Hertzian contact force arises only when there is contact between the ball and the raceway, the respective contact force is set to zero when the contact deformation is equal or smaller than zero. This is indicated by subscript "+". θ i is given by where ω c is the speed of the cage.
where ω is the shaft speed. Based on the structural dynamics theory, the equations of motion for a two-degrees-of-freedom system can be written as follows: m .. where m is the mass of the rotor and of the inner race, c is the damping factor, and W is the static load. The damping factor c is commonly chosen to be 200 Ns/m [24,25,[27][28][29].

Load Dependent Friction Torque
A reasonable estimate of the friction torque of a given rolling bearing under moderate load and speed conditions is the sum of load friction torque M 1 and viscous friction torque M v . M v is primarily dependent on speed and the method of lubrication [23] and is independent of bearing load. Hence, M v isn't going to be discussed because of its insensitivity to bearing defects in this paper.
The load dependent friction torque is given by where f 1 is a factor dependent on the bearing design, and F β depends on the magnitude and direction of the applied load. The detailed formulas and parameters for f 1 are introduced by literature [23]. For ball bearings having a nominal contact angle 0 • , there is not applied axial load, so F β equals the applied radial load F r . For static equilibrium, the applied radial load must equal the sum of the Hertzian forces from each of the balls.

Simulated Hertzian Forces and Friction Torque
The spall model and dynamic model of the bearing were solved in MATLAB using the ordinary differential equation solver (ode45). For the defective SKF 6206 bearing in a 2.2 kW induction motor, the parameters of the spall model and dynamic model of the bearing are given in Table 1.  Table 1, the number of balls N, ball radius r b , pitch radius r c , radial clearance C r and load-deflection factor K are parameters related to the SKF 6206 bearing. The defect depth h d , angle of the defect centre Φ d and angle of the defect length ∆Φ d are parameters related to the assumed spall in bearing outer raceway. The mass of the rotor m and static load W are parameters related to the 2.2 kW induction motor, and the W is equal to m times gravity acceleration. For the bearing SKF 6206 mounted in the motor, the damping factor c is commonly chosen to be 200 Ns/m. The motor speed n and shaft speed ω are parameters related to the motor operation, and ω = 2πn/60.
The total Hertzian forces along the X-axis and Y-axis are as shown in Figure 5. In Figure 5, because Φd is 0°, which means that the spall is just located in the 6 o'clock region, the X-axis total Hertzian forces FX is much larger than the FY, which will be ignored in following discussion.
Whether in FX or in FY, the cycle of the Hertzian forces caused by the spall exactly corresponds to the characteristic vibration frequencies fc, which are given by: where fr is the mechanical rotor frequency, N is the number of balls, Db is the diameter of the ball, Dc is bearing pitch diameter, and α is the contact angle between the ball and the raceway.
Except for fc, the information about spall size is also revealed distinctly in Figure 5. There are two impulse responses, respectively located at the entry and exit of the spall. The spall size corresponds to the primary interval between "d" and "e", the additional interval between "b" and "d", and the additional interval between "e" and "g". The primary interval is easy to determine based on Figure 5. The latter two additional intervals are not very distinct, but can be estimated by size of the ball.
Through different pathways, two impulse responses caused by the entry and exit of the spall can also be propagated to the vibration signal and the acoustic emission signal. Then, spall size estimation can be implemented by determining the average interval between two events [2,6,[9][10][11][12][13][14].
Whereafter, the ΔΦd is set as 5°, 10° and 15° respectively in Table 1; the friction torques are as shown in Figure 6.
In Figure 6, the two impulse responses caused by the two events exhibit consistent in friction torque for different spall lengths, but the amplitude of the friction torque is very small. For spall size estimation based on current, not only the friction torque will be buried, likely in load torque, but also the friction torque modulated in current will be faced with strong noise which will include the power supply frequency and its harmonic. All this means that spall size estimation based on current is scarcely possible by determining the average interval in the time series. As an alternative, in the frequency domain, tracking the frequency feature of the friction torque modulated in current will be a feasible solution. In Figure 5, because Φ d is 0 • , which means that the spall is just located in the 6 o'clock region, the X-axis total Hertzian forces F X is much larger than the F Y , which will be ignored in following discussion.
Whether in F X or in F Y , the cycle of the Hertzian forces caused by the spall exactly corresponds to the characteristic vibration frequencies f c , which are given by: where f r is the mechanical rotor frequency, N is the number of balls, D b is the diameter of the ball, D c is bearing pitch diameter, and α is the contact angle between the ball and the raceway. Except for f c , the information about spall size is also revealed distinctly in Figure 5. There are two impulse responses, respectively located at the entry and exit of the spall. The spall size corresponds to the primary interval between "d" and "e", the additional interval between "b" and "d", and the additional interval between "e" and "g". The primary interval is easy to determine based on Figure 5. The latter two additional intervals are not very distinct, but can be estimated by size of the ball.
Through different pathways, two impulse responses caused by the entry and exit of the spall can also be propagated to the vibration signal and the acoustic emission signal. Then, spall size estimation can be implemented by determining the average interval between two events [2,6,[9][10][11][12][13][14].
Whereafter, the ∆Φ d is set as 5 • , 10 • and 15 • respectively in Table 1; the friction torques are as shown in Figure 6.
In Figure 6, the two impulse responses caused by the two events exhibit consistent in friction torque for different spall lengths, but the amplitude of the friction torque is very small. For spall size estimation based on current, not only the friction torque will be buried, likely in load torque, but also the friction torque modulated in current will be faced with strong noise which will include the power supply frequency and its harmonic. All this means that spall size estimation based on current is scarcely possible by determining the average interval in the time series. As an alternative, in the frequency domain, tracking the frequency feature of the friction torque modulated in current will be a feasible solution.
Unlike a single-point defect, which only produces a single characteristic vibration frequency, the spall including two impulse responses will produce multiple characteristic vibration frequencies. Therefore, the relation between spall length and the multiple characteristic vibration frequencies in friction torque must be investigated first. Unlike a single-point defect, which only produces a single characteristic vibration frequency, the spall including two impulse responses will produce multiple characteristic vibration frequencies. Therefore, the relation between spall length and the multiple characteristic vibration frequencies in friction torque must be investigated first.

The Relation Between Spall Size and Multiple Characteristic Vibration Frequencies
For constant spall depth hd = 100 μm, Changing the ΔΦd from 1° to 20° with a distance of 1°, the friction torque was respectively solved. For each friction torque, the amplitudes of 1st~3rd multiple characteristic vibration frequencies fc, 2fc and 3fc were achieved based on Fast Fourier Transform (FFT) in MATLAB. Only the fc, 2fc and 3fc frequencies are analyzed in this paper because of the sufficient spall length information carried by these three frequencies, even though other, higher, multiplecharacteristic vibration frequencies are also present.
When the motor speed was set as 150, 300 and 450 r/min respectively, the amplitude trend of fc, 2fc and 3fc in friction torque for various spall length were as shown in Figure 7.

The Relation Between Spall Size and Multiple Characteristic Vibration Frequencies
For constant spall depth h d = 100 µm, Changing the ∆Φ d from 1 • to 20 • with a distance of 1 • , the friction torque was respectively solved. For each friction torque, the amplitudes of 1st~3rd multiple characteristic vibration frequencies f c , 2f c and 3f c were achieved based on Fast Fourier Transform (FFT) in MATLAB. Only the f c , 2f c and 3f c frequencies are analyzed in this paper because of the sufficient spall length information carried by these three frequencies, even though other, higher, multiple-characteristic vibration frequencies are also present.
When the motor speed was set as 150, 300 and 450 r/min respectively, the amplitude trend of f c , 2f c and 3f c in friction torque for various spall length were as shown in Figure 7. Unlike a single-point defect, which only produces a single characteristic vibration frequency, the spall including two impulse responses will produce multiple characteristic vibration frequencies. Therefore, the relation between spall length and the multiple characteristic vibration frequencies in friction torque must be investigated first.

The Relation Between Spall Size and Multiple Characteristic Vibration Frequencies
For constant spall depth hd = 100 μm, Changing the ΔΦd from 1° to 20° with a distance of 1°, the friction torque was respectively solved. For each friction torque, the amplitudes of 1st~3rd multiple characteristic vibration frequencies fc, 2fc and 3fc were achieved based on Fast Fourier Transform (FFT) in MATLAB. Only the fc, 2fc and 3fc frequencies are analyzed in this paper because of the sufficient spall length information carried by these three frequencies, even though other, higher, multiplecharacteristic vibration frequencies are also present.
When the motor speed was set as 150, 300 and 450 r/min respectively, the amplitude trend of fc, 2fc and 3fc in friction torque for various spall length were as shown in Figure 7.  In Figure 7, following the increase of the spall length, the fluctuation of amplitude of 2fc and 3fc are more obvious than fc. When ΔΦd reached a certain level (10° when motor speed is 300r/min shown in Figure 7b), the amplitude of 3fc reached its maximum. When the ΔΦd further increased (15°-16°), the amplitude of 3fc reduced to the minimum; in the meantime, the amplitude of 2fc reached its maximum. The regular rise and fall of amplitude of 2fc and 3fc was consistent for different motor speeds. When the motor speed is slower, the amplitude trend collectively shifts left, as shown in Figure 7a. When the motor speed is higher, the amplitude trend collectively shifts right, as shown in Figure 7c.
Following an increase in the motor speed, the amplitudes of fc, 2fc and 3fc all increase. Theoretically, increased motor speed is possibly beneficial for identifying those spall-related characteristic frequencies in stator current.
In order to verify that the regular rise and fall of amplitude of 2fc and 3fc is universal for other cases, the various spall depth is also discussed. When the motor speed was set as 300 r/min, for three spall depths hd = 20, 40 and 60 μm respectively, the amplitude trends of fc, 2fc and 3fc were as shown in Figure 8.    Figure 7, following the increase of the spall length, the fluctuation of amplitude of 2f c and 3f c are more obvious than f c . When ∆Φ d reached a certain level (10 • when motor speed is 300r/min shown in Figure 7b), the amplitude of 3f c reached its maximum. When the ∆Φ d further increased (15 • -16 • ), the amplitude of 3f c reduced to the minimum; in the meantime, the amplitude of 2f c reached its maximum. The regular rise and fall of amplitude of 2f c and 3f c was consistent for different motor speeds. When the motor speed is slower, the amplitude trend collectively shifts left, as shown in Figure 7a. When the motor speed is higher, the amplitude trend collectively shifts right, as shown in Figure 7c.
Following an increase in the motor speed, the amplitudes of f c , 2f c and 3f c all increase. Theoretically, increased motor speed is possibly beneficial for identifying those spall-related characteristic frequencies in stator current.
In order to verify that the regular rise and fall of amplitude of 2f c and 3f c is universal for other cases, the various spall depth is also discussed. When the motor speed was set as 300 r/min, for three spall depths h d = 20, 40 and 60 µm respectively, the amplitude trends of f c , 2f c and 3f c were as shown in Figure 8. In Figure 7, following the increase of the spall length, the fluctuation of amplitude of 2fc and 3fc are more obvious than fc. When ΔΦd reached a certain level (10° when motor speed is 300r/min shown in Figure 7b), the amplitude of 3fc reached its maximum. When the ΔΦd further increased (15°-16°), the amplitude of 3fc reduced to the minimum; in the meantime, the amplitude of 2fc reached its maximum. The regular rise and fall of amplitude of 2fc and 3fc was consistent for different motor speeds. When the motor speed is slower, the amplitude trend collectively shifts left, as shown in Figure 7a. When the motor speed is higher, the amplitude trend collectively shifts right, as shown in Figure 7c.
Following an increase in the motor speed, the amplitudes of fc, 2fc and 3fc all increase. Theoretically, increased motor speed is possibly beneficial for identifying those spall-related characteristic frequencies in stator current.
In order to verify that the regular rise and fall of amplitude of 2fc and 3fc is universal for other cases, the various spall depth is also discussed. When the motor speed was set as 300 r/min, for three spall depths hd = 20, 40 and 60 μm respectively, the amplitude trends of fc, 2fc and 3fc were as shown in Figure 8.  In Figure 8, for different spall depths, the regular rise and fall of amplitude of 2f c and 3f c is still subsistent. In addition, following an increase in the spall depth, the amplitudes of f c , 2f c and 3f c just increase a little. That is to say, the amplitude of the multiple characteristic vibration frequencies is not very sensitive for spall depth.
Based on termly observation of current during the motor operation, the spall length may be estimated approximately by tracking the regular rise and fall of amplitude of 2f c and 3f c . The estimation of spall length will be sufficient evidence for the correct maintenance decision and a reasonable maintenance schedule for the induction motor.

The Multiple Characteristic Vibration Frequencies Modulated in Current
The friction torque is transmitted to the stator current by means of nonlinear phase modulation, which, even though a single frequency, can also be evolved into a chain of frequencies in current. If there are spectrum overlaps among the multiple characteristic vibration frequencies modulated in the current, spall size estimation based on the regular rise and fall of amplitude of 2f c and 3f c will be imprecise. Therefore, the multiple characteristic vibration frequencies modulated in current will be investigated.
Only taking account of the effect of f c , 2f c and 3f c , which are the primary parts of fluctuant friction torque caused by spall, the total load torque T l can be described as follows: where T 0 is the constant torque, n is a positive integer, ω c is the characteristic vibration angular frequency and T n is the amplitude of nω c . Based on the same method introduced in [36], the stator current is given by where ω s is the electrical supply angular frequency, β n is the modulation index of the n th multiple characteristic vibration angular frequency nω c , and β n is given by where p is the pole pairs and J is the total inertia of the rotating system. Following an increase of n, β n will decrease. After modulating in current, the amplitude of 2f c will be less than that shown in Figure 7; the amplitude of 3f c will also be lesser, but the regular rise and fall is unaltered. For this reason, those other higher order multiple characteristic vibration frequencies are not taken into account, because those modulated in current are scarcely possible to identify.
Equation (17) can be expanded on the basis of the first-kind Bessel function where the symbol "Im" represents taking the imaginary part, m is an integer, and J m (β n ) is the m th Bessel function for β n .
For n=1, 2, and 3, Equation (19) can be simplified as where a, b, and c are integers, and similar to m in (19). J a (β 1 ) is the a th Bessel function for β 1 , J b (β 2 ) is the b th Bessel function for β 2 , and J c (β 3 ) is the c th Bessel function for β 3 .
Letting k = a + 2b + 3c, then, the bearing fault-related frequency f bf in stator current is given by where f s is the power supply frequency; f s ± kf c is commonly called the k th sideband in phase modulation. Equation (21) is same as the proposed result in the majority of literature, but the difference is that the amplitude of side frequencies in current can be revealed using Equation (20).
For the general formula of J a (β 1 ), J b (β 2 ), and J c (β 3 ), J m (β n ) can be expressed as Considering that equation (22) can be rewritten as Considering that the magnitude of β n is extremely small (less than 0.1) [3,16,36], only k = 0 deserves attention. This therefore leads to Equation (25) can be directly rewritten as In addition, the symmetry of the Bessel function is shown, as J -m (β n ) = −J m (β n ) when m is odd, and J -m (β n ) = J m (β n ) when m is even.
Focusing on the amplitude of the first sideband f s +f c , k=1, that is, a+2b+3c=1. There are innumerable solutions for a, b, and c.
Considering that β n is extremely small, the magnitudes of 0.25I 1 β 1 β 2 and 0.25I 1 β 2 β 3 are much smaller than the magnitude of 0.5I 1 β 1, and are ignorable. The amplitude of f s + f c is close to 0.5I 1 β 1 .
In the same way, the amplitude of other sidebands can be calculated. The final conclusion is that the amplitude of f s ± f c is 0.5I 1 β 1 ; the amplitude of f s ± 2f c is 0.5I 1 β 2 ; the amplitude of f s ± 3f c is 0.5I 1 β 3 .
This conclusion means that, after transmitting to current, the spall size information implied in the 1st~3rd multiple characteristic vibration frequencies presents perfectly in the current without spectrum overlap. In brief, the amplitude of f s ± 2f c is overwhelmingly contributed by the 2nd multiple characteristic vibration frequencies, 2f c , and is unaffected by the other multiple characteristic vibration frequencies.
The condition of f s ± f c and f s ± 3f c is similar to that of f s ± 2f c .
In theory, spall size estimation based on current can be implemented by tracking the regular rise and fall of amplitude of multiple characteristic vibration frequencies modulated in current.

Experimental Verification on Fluctuant Friction Torque
Based on the simulation results in Section 2.5, there are two impulse responses respectively located at the entry and exit of the spall, which present in friction torque. To verify this simulation result, experimental verification is necessary. However, in practice, the friction torque of a spalled bearing assembled in a motor, is hard to be measured directly. Therefore, an indirect experimental scheme is implemented for only the spalled bearing instead of the bearing assembled in motor. Theoretically, because there are differences in just two parameters, mass of the rotor and static load, the amplitudes of two impulse responses will be weakened overall, but the forms of two impulse responses stay the same.
The M9908B bearing friction torque measuring instrument produced by Luoyang Bearing Science & Technology Co. Ltd in China, was used in the present research. The M9908B is as shown in Figure 9. This conclusion means that, after transmitting to current, the spall size information implied in the 1st~3rd multiple characteristic vibration frequencies presents perfectly in the current without spectrum overlap. In brief, the amplitude of fs ± 2fc is overwhelmingly contributed by the 2nd multiple characteristic vibration frequencies, 2fc, and is unaffected by the other multiple characteristic vibration frequencies. The condition of fs ± fc and fs ± 3fc is similar to that of fs ± 2fc.
In theory, spall size estimation based on current can be implemented by tracking the regular rise and fall of amplitude of multiple characteristic vibration frequencies modulated in current.

Experimental Verification on Fluctuant Friction Torque
Based on the simulation results in Section 2.5, there are two impulse responses respectively located at the entry and exit of the spall, which present in friction torque. To verify this simulation result, experimental verification is necessary. However, in practice, the friction torque of a spalled bearing assembled in a motor, is hard to be measured directly. Therefore, an indirect experimental scheme is implemented for only the spalled bearing instead of the bearing assembled in motor. Theoretically, because there are differences in just two parameters, mass of the rotor and static load, the amplitudes of two impulse responses will be weakened overall, but the forms of two impulse responses stay the same.
The M9908B bearing friction torque measuring instrument produced by Luoyang Bearing Science & Technology Co. Ltd in China, was used in the present research. The M9908B is as shown in Figure 9. The technical specifications of the M9908B are as follows: torque range is 50 mN.m; rotational speed of drive spindle is 5 rpm; error is ± 0.5 mN.m; repeatability is ± 0.5 mN.m.
Like the majority of bearing friction torque measuring instruments (for example Vibrac BRG 3000), an axial load is commonly used in the testing of bearings with M9908B. It is clear that the fluctuation of friction torque caused by spall cannot be perfectly shown by means of axial loading. Hence, a customized radial loading unit was supplemented to M9908B.
The bearing's outer ring is driven by the additional radial loading unit fixed in the drive spindle. The bearing's inner ring is fixed in the air spindle, connected with a torque transducer. When the bearing's outer ring rotates along with the drive spindle, because of the presence of the friction torque, the bearing's inner ring shows a tendency of rotation along with the outer ring. The rotation tendency of the inner ring will be impeded by a torque transducer; for dynamic equilibrium, the friction torque was measured. The technical specifications of the M9908B are as follows: torque range is 50 mN.m; rotational speed of drive spindle is 5 rpm; error is ± 0.5 mN.m; repeatability is ± 0.5 mN.m.
Like the majority of bearing friction torque measuring instruments (for example Vibrac BRG 3000), an axial load is commonly used in the testing of bearings with M9908B. It is clear that the fluctuation of friction torque caused by spall cannot be perfectly shown by means of axial loading. Hence, a customized radial loading unit was supplemented to M9908B.
The bearing's outer ring is driven by the additional radial loading unit fixed in the drive spindle. The bearing's inner ring is fixed in the air spindle, connected with a torque transducer. When the bearing's outer ring rotates along with the drive spindle, because of the presence of the friction torque, the bearing's inner ring shows a tendency of rotation along with the outer ring. The rotation tendency of the inner ring will be impeded by a torque transducer; for dynamic equilibrium, the friction torque was measured.
Basing on electro-discharge machining, the outer raceway of the SKF 6206 bearings was artificially drilled with an approximate angle of spall length of 5 • , 10 • and 15 • respectively, as shown in Figure 10. Basing on electro-discharge machining, the outer raceway of the SKF 6206 bearings was artificially drilled with an approximate angle of spall length of 5°, 10° and 15° respectively, as shown in Figure 10. For the above three defective bearings, the friction torques achieved respectively based on M9908B were as shown in Figure 11. As Figure 11 shows, the simulation results in Section 2.5 were verified experimentally: 1. The cycle of the fluctuant friction torque exactly corresponded to the characteristic vibration frequencies fc. 2. There were two impulse responses when the ball passed through the spall. 3. The interval between the two impulse responses approximatively corresponded to the spall length.
Therefore, the regular rise and fall of amplitude of 2fc and 3fc in friction torque caused by the interval between the two impulse responses is feasible for estimating spall length.

Experimental Platform
For the above three defective bearings, the friction torques achieved respectively based on M9908B were as shown in Figure 11. Basing on electro-discharge machining, the outer raceway of the SKF 6206 bearings was artificially drilled with an approximate angle of spall length of 5°, 10° and 15° respectively, as shown in Figure 10. For the above three defective bearings, the friction torques achieved respectively based on M9908B were as shown in Figure 11. As Figure 11 shows, the simulation results in Section 2.5 were verified experimentally: 1. The cycle of the fluctuant friction torque exactly corresponded to the characteristic vibration frequencies fc. 2. There were two impulse responses when the ball passed through the spall. 3. The interval between the two impulse responses approximatively corresponded to the spall length.
Therefore, the regular rise and fall of amplitude of 2fc and 3fc in friction torque caused by the interval between the two impulse responses is feasible for estimating spall length. As Figure 11 shows, the simulation results in Section 2.5 were verified experimentally:
The cycle of the fluctuant friction torque exactly corresponded to the characteristic vibration frequencies f c .

2.
There were two impulse responses when the ball passed through the spall. 3.
The interval between the two impulse responses approximatively corresponded to the spall length.
Therefore, the regular rise and fall of amplitude of 2f c and 3f c in friction torque caused by the interval between the two impulse responses is feasible for estimating spall length.

Experimental Platform
An experimental platform was designed based on FPGA and IPM. The schematic diagram is shown as Figure 12. An experimental platform was designed based on FPGA and IPM. The schematic diagram is shown as Figure 12. In Figure 12, T is a 3-phase autotransformer; R is a 3-phase rectifier; IPM is a Mitsubishi Intelligent Power Modules PM50CSD120, with a collector current of 50 A and collector-emitter voltage of 1200V; the Optic coupler is a HCPL 4504; the drive circuit is a triode 9013; FPGA is an Altera Cyclone IV E devices EP4CE15F17C8. SPWM modulation mode was employed in the FPGA drive system operating at carrier frequency 20 kHz.
The parameters of three-phase induction motor were as follows: rated power 2.2 kW; rated voltage 380 V; rated frequency 50 Hz; rated current 5.03 A; rated speed 1430 r/min; pole pairs 2; bearing type SKF 6206.
Induction motor speed was measured by handheld tachometer. Induction motor stator current was measured by Hall current sensors LEM IT60-S, whose measuring range was 60A. Linearity error was 20 ppm, and frequency bandwidth was 500 kHz.
The current data was sampled by National Instruments X series multifunction data acquisition card PXIe-6356, with an analog digital converter resolution of 16 bits, sampling rate of 1.25 MS/s， and bandwidth of 1 MHz. PC is a National Instruments 2.3 GHz quad-core embedded controller PXIe-8135. A picture of the experimental platform is shown as Figure 13.
Theoretically, the lower the motor speed is, the more distinct the two impulse responses in friction torque are, without overlap between two impulse responses. However, in practice, when the power supply frequency is set too low, the fault-related frequencies probably overlap with power supply frequency and its harmonic, because of the wide lobe produced by spectrum analysis. Therefore, the power supply frequency is initially set to 10 Hz.
In order to avoid false positive detection caused by the fluctuation of mechanical load or shaft coupling, the motor was in no-load.
To ensure reasonable comparison among spectrum analysis results, the Welch method was adopted for spectrum analysis throughout the following research. The specific parameters were as follows: FFT points were 8192, the sampling points were 8192, sampling frequency was 1000 Hz. In Figure 12, T is a 3-phase autotransformer; R is a 3-phase rectifier; IPM is a Mitsubishi Intelligent Power Modules PM50CSD120, with a collector current of 50 A and collector-emitter voltage of 1200V; the Optic coupler is a HCPL 4504; the drive circuit is a triode 9013; FPGA is an Altera Cyclone IV E devices EP4CE15F17C8. SPWM modulation mode was employed in the FPGA drive system operating at carrier frequency 20 kHz.
The parameters of three-phase induction motor were as follows: rated power 2.2 kW; rated voltage 380 V; rated frequency 50 Hz; rated current 5.03 A; rated speed 1430 r/min; pole pairs 2; bearing type SKF 6206.
Induction motor speed was measured by handheld tachometer. Induction motor stator current was measured by Hall current sensors LEM IT60-S, whose measuring range was 60A. Linearity error was 20 ppm, and frequency bandwidth was 500 kHz.
The current data was sampled by National Instruments X series multifunction data acquisition card PXIe-6356, with an analog digital converter resolution of 16 bits, sampling rate of 1.25 MS/s, and bandwidth of 1 MHz. PC is a National Instruments 2.3 GHz quad-core embedded controller PXIe-8135. A picture of the experimental platform is shown as Figure 13.

Highlighting Fault Features Based on Reduced Voltage Frequency Ratio
For tracking multiple characteristic vibration frequencies modulated in current, the squared envelope of current based on Hilbert transform is employed to suppress noise. However, this method alone is far from enough; therefore, this paper attempts to improve the SNR based on additional approaches.
To supply enough electromagnetic torque to the load, the motor commonly operates in rated Theoretically, the lower the motor speed is, the more distinct the two impulse responses in friction torque are, without overlap between two impulse responses. However, in practice, when the power supply frequency is set too low, the fault-related frequencies probably overlap with power supply frequency and its harmonic, because of the wide lobe produced by spectrum analysis. Therefore, the power supply frequency is initially set to 10 Hz.
In order to avoid false positive detection caused by the fluctuation of mechanical load or shaft coupling, the motor was in no-load.
To ensure reasonable comparison among spectrum analysis results, the Welch method was adopted for spectrum analysis throughout the following research. The specific parameters were as follows: FFT points were 8192, the sampling points were 8192, sampling frequency was 1000 Hz.

Highlighting Fault Features Based on Reduced Voltage Frequency Ratio
For tracking multiple characteristic vibration frequencies modulated in current, the squared envelope of current based on Hilbert transform is employed to suppress noise. However, this method alone is far from enough; therefore, this paper attempts to improve the SNR based on additional approaches.
To supply enough electromagnetic torque to the load, the motor commonly operates in rated main flux, in which the voltage frequency ratio is always constant in common inverter. Because the amplitude of friction torque is much less than the amplitude of electromagnetic torque, the fluctuant friction torque is filtered rapidly by motor inertia before modulating in current, which dynamic process is very transient. Accordingly, the fault features modulated in current are very weak.
Under the constant power supply frequency without any form of closed-loop control, the main flux and the electromagnetic torque of the motor will be reduced following the reduction of supply voltage; meanwhile, the fluctuant friction torque caused by spall is almost unchanged. When the supply voltage drops to a certain extent, the electromagnetic torque is also close to the friction torque. In which case, when the ball passes through the spall, the fault features will be preferably modulated in the current due to the longish dynamic interactive process between the fluctuant friction torque and the electromagnetic torque, together with the slower motor speed. Consequently, the fault features modulated in the current will be highlighted.
In the experimental platform in Figure 13, the power supply frequency was held constant at 10 Hz, the motor was in no-load, and the bearing had a spall length of 5 • (Figure 10a). For the experimental motor, normal voltage frequency ratio was 380 V/50 Hz (7.6). The voltage frequency ratio was reduced gradually for verifying the effect of changed main flux based on the FPGA and IPM drive system. During experiment, the voltage frequency ratio was set as 76 V/10 Hz (7.6), 38 V/10 Hz (3.8), 19 V/10 Hz (1.9) respectively; the measured motor speed was 299 r/min, 298 r/min, 293 r/min, and the f c was 17.7 Hz, 17.7 Hz, 17.4 Hz accordingly. The squared envelope spectra of the stator currents in all of above cases were as shown as Figure 14.

Highlighting Fault Features Based on Reduced Voltage Frequency Ratio
For tracking multiple characteristic vibration frequencies modulated in current, the squared envelope of current based on Hilbert transform is employed to suppress noise. However, this method alone is far from enough; therefore, this paper attempts to improve the SNR based on additional approaches.
To supply enough electromagnetic torque to the load, the motor commonly operates in rated main flux, in which the voltage frequency ratio is always constant in common inverter. Because the amplitude of friction torque is much less than the amplitude of electromagnetic torque, the fluctuant friction torque is filtered rapidly by motor inertia before modulating in current, which dynamic process is very transient. Accordingly, the fault features modulated in current are very weak.
Under the constant power supply frequency without any form of closed-loop control, the main flux and the electromagnetic torque of the motor will be reduced following the reduction of supply voltage; meanwhile, the fluctuant friction torque caused by spall is almost unchanged. When the supply voltage drops to a certain extent, the electromagnetic torque is also close to the friction torque. In which case, when the ball passes through the spall, the fault features will be preferably modulated in the current due to the longish dynamic interactive process between the fluctuant friction torque and the electromagnetic torque, together with the slower motor speed. Consequently, the fault features modulated in the current will be highlighted.
In the experimental platform in Figure 13, the power supply frequency was held constant at 10 Hz, the motor was in no-load, and the bearing had a spall length of 5° (Figure 10a). For the experimental motor, normal voltage frequency ratio was 380 V/50 Hz (7.6). The voltage frequency ratio was reduced gradually for verifying the effect of changed main flux based on the FPGA and IPM drive system. During experiment, the voltage frequency ratio was set as 76 V/10 Hz (7.6), 38 V/10 Hz (3.8), 19 V/10 Hz (1.9) respectively; the measured motor speed was 299 r/min, 298 r/min, 293 r/min, and the fc was 17.7 Hz, 17.7 Hz, 17.4 . In Figure 14a, when the voltage frequency ratio was normal (7.6), there were four fault frequencies: 2fs + fc, 6fs − fc, 4fs + fc, 8fs − fc. All these frequencies were below the average of spectral lines, and therefore are not suitable for use as a diagnosis index. A similar situation is also presented in Figure 14b, where the voltage frequency ratio was half of normal (3.8).
Following the reduction of voltage frequency ratio further, the motor entered the critical state between rotation and non-rotation. For the experimental motor, when the voltage frequency ratio was a quarter of normal (1.9), the motor just held a steady rotation speed. In Figure 14c, the four fault frequencies: 2fs + fc, 6fs − fc, 4fs + fc, 8fs − fc are unusually noticeable. Consequently, the objective of enhancing the SNR of defect signatures is achieved perfectly.
In order to verify the performance of proposed method, some traditional approaches were also used to compare. The instantaneous symmetrical components and Park's vector are often used to highlight fault features based on three phase currents. When the voltage frequency ratio was normal, the squared envelope spectra based on traditional approaches were as shown as Figure 15. In Figure 15a, there are four fault frequencies: fs + fc, 5fs − fc, 6fs − fc and 4fs + fc are identifiable. In Figure 15b, there are three fault frequencies: fs + fc, 6fs − fc and 4fs + fc are identifiable. All these frequencies are more distinct than those in Figure 14a. Therefore, the fault features may be highlighted based on instantaneous symmetrical components and Park's vector.
However, the performance of fault feature highlighting based on traditional approaches is far inferior to the proposed method, as is shown in Figure 14c. In fact, because the motor is in no-load, the defect signatures are inappreciable under low bearing radial load, which was always a challenging issue facing low load up to now. Nevertheless, the proposed method addresses this issue In Figure 14a, when the voltage frequency ratio was normal (7.6), there were four fault frequencies: 2f s + f c , 6f s − f c , 4f s + f c , 8f s − f c . All these frequencies were below the average of spectral lines, and therefore are not suitable for use as a diagnosis index. A similar situation is also presented in Figure 14b, where the voltage frequency ratio was half of normal (3.8).
Following the reduction of voltage frequency ratio further, the motor entered the critical state between rotation and non-rotation. For the experimental motor, when the voltage frequency ratio was a quarter of normal (1.9), the motor just held a steady rotation speed. In Figure 14c, the four fault frequencies: 2f s + f c , 6f s − f c , 4f s + f c , 8f s − f c are unusually noticeable. Consequently, the objective of enhancing the SNR of defect signatures is achieved perfectly.
In order to verify the performance of proposed method, some traditional approaches were also used to compare. The instantaneous symmetrical components and Park's vector are often used to highlight fault features based on three phase currents. When the voltage frequency ratio was normal, the squared envelope spectra based on traditional approaches were as shown as Figure 15. In Figure 14a, when the voltage frequency ratio was normal (7.6), there were four fault frequencies: 2fs + fc, 6fs − fc, 4fs + fc, 8fs − fc. All these frequencies were below the average of spectral lines, and therefore are not suitable for use as a diagnosis index. A similar situation is also presented in Figure 14b, where the voltage frequency ratio was half of normal (3.8).
Following the reduction of voltage frequency ratio further, the motor entered the critical state between rotation and non-rotation. For the experimental motor, when the voltage frequency ratio was a quarter of normal (1.9), the motor just held a steady rotation speed. In Figure 14c, the four fault frequencies: 2fs + fc, 6fs − fc, 4fs + fc, 8fs − fc are unusually noticeable. Consequently, the objective of enhancing the SNR of defect signatures is achieved perfectly.
In order to verify the performance of proposed method, some traditional approaches were also used to compare. The instantaneous symmetrical components and Park's vector are often used to highlight fault features based on three phase currents. When the voltage frequency ratio was normal, the squared envelope spectra based on traditional approaches were as shown as Figure 15. In Figure 15a, there are four fault frequencies: fs + fc, 5fs − fc, 6fs − fc and 4fs + fc are identifiable. In Figure 15b, there are three fault frequencies: fs + fc, 6fs − fc and 4fs + fc are identifiable. All these frequencies are more distinct than those in Figure 14a. Therefore, the fault features may be highlighted based on instantaneous symmetrical components and Park's vector.
However, the performance of fault feature highlighting based on traditional approaches is far inferior to the proposed method, as is shown in Figure 14c. In fact, because the motor is in no-load, the defect signatures are inappreciable under low bearing radial load, which was always a challenging issue facing low load up to now. Nevertheless, the proposed method addresses this issue In Figure 15a, there are four fault frequencies: f s + f c , 5f s − f c , 6f s − f c and 4f s + f c are identifiable. In Figure 15b, there are three fault frequencies: f s + f c , 6f s − f c and 4f s + f c are identifiable. All these frequencies are more distinct than those in Figure 14a. Therefore, the fault features may be highlighted based on instantaneous symmetrical components and Park's vector.
However, the performance of fault feature highlighting based on traditional approaches is far inferior to the proposed method, as is shown in Figure 14c. In fact, because the motor is in no-load, the defect signatures are inappreciable under low bearing radial load, which was always a challenging issue facing low load up to now. Nevertheless, the proposed method addresses this issue by reducing the electromagnetic torque to a certain extent, just matching weak friction torque. In this case, it is inevitable that the fault features will be highlighted perfectly.
Furthermore, the proposed method may be also applicable to loaded motors without any form of closed-loop control.

The Effect of Rotation Speed on Current Sideband
The amplitudes of characteristic vibration frequencies increase slightly following an increased motor speed, based on the theoretical analysis in Section 2.6. The effect of rotation speed on current sideband will be investigated in the present section.
On one hand, the optimal power supply frequency must be selected. On the other hand, the multiple characteristic vibration frequencies need to be presented collectively for intercomparison. Considering the above two points, the bearing with spall length of 10 • (Figure 10b), is employed in experimental motor, the motor is in no-load, and the voltage frequency ratio is set as a quarter of normal (1.9) throughout.
During the experiment, the power supply frequencies were set as 6 Hz, 10 Hz, 12 Hz, 15 Hz and 20 Hz respectively; the measured motor speeds were 162 r/min, 290 r/min, 354 r/min, 444 r/min and 589 r/min, and f c values were 9.6 Hz, 17.2 Hz, 21 Hz, 26.3 Hz and 34.9 Hz accordingly. The squared envelope spectra of stator currents in all of above cases were as shown as Figure 16.
Sensors 2020, 20, 1631 18 of 23 by reducing the electromagnetic torque to a certain extent, just matching weak friction torque. In this case, it is inevitable that the fault features will be highlighted perfectly. Furthermore, the proposed method may be also applicable to loaded motors without any form of closed-loop control.

The Effect of Rotation Speed on Current Sideband
The amplitudes of characteristic vibration frequencies increase slightly following an increased motor speed, based on the theoretical analysis in Section 2.6. The effect of rotation speed on current sideband will be investigated in the present section.
On one hand, the optimal power supply frequency must be selected. On the other hand, the multiple characteristic vibration frequencies need to be presented collectively for intercomparison. Considering the above two points, the bearing with spall length of 10° (Figure 10b), is employed in experimental motor, the motor is in no-load, and the voltage frequency ratio is set as a quarter of normal (1.9) throughout.
During the experiment, the power supply frequencies were set as 6 Hz, 10 Hz, 12 Hz, 15 Hz and 20 Hz respectively; the measured motor speeds were 162 r/min, 290 r/min, 354 r/min, 444 r/min and 589 r/min, and fc values were 9.6 Hz, 17.2 Hz, 21 Hz, 26.3 Hz and 34.9 Hz accordingly. The squared envelope spectra of stator currents in all of above cases were as shown as Figure 16.  In Figure 16a, the power supply frequency is 6 Hz, a number of fault frequencies arise, far beyond the average of spectral lines. The sidebands of 3f c included: 6f s +3f c , 4f s +3f c and f s +3f c . The sidebands of 2f c included: 4f s +2f c and 3f s +2f c . The sidebands of f c included: f s +f c , 3f s +f c , 2f s +f c , 3f s −f c and 4f s +f c . Among all fault frequencies, 6f s +3f c , modulated from 3f c , was maximal, next was f s +f c and 4f s +2f c . Theoretically, this is a suitable power supply frequency for investigating the multiple modulated frequencies. However, the modulated frequencies were too close to each other, and several modulated frequencies even overlapped with the power supply frequency and its harmonics.
In Figure 16b, the power supply frequency was 10Hz, and number of fault frequencies arose similarly. For simplicity, no more enumerating every frequency; only the most prominent frequency among each group frequencies will be chosen for comparison. Obviously, 4f s −3f c was maximal, next was 5f s −f c and 2f s +2f c . Because there was no sideband overlap with the power supply frequency and its harmonics, this is a very suitable power supply frequency for investigating the multiple modulated frequencies.
In Figure 16c, the power supply frequency was 12 Hz. f s +3f c was maximal, next was 2f s + f c and 7f s −2f c . On the whole, the amplitudes of all fault frequencies were slightly lower than above two cases.
In Figure 16d, the power supply frequency was 15 Hz. f s +3f c was maximal, next was 3f s +f c and 5f s −2f c . For the sidebands of 2f c , except the prominent frequency 5f s −2f c , there was only one frequency, f s −2f c , which is almost masked in spectral lines. Some fault frequencies were also not prominent, as in the above three cases.
In Figure 16e, the power supply frequency was 20Hz. 2f s +f c was maximal, next was 2f s −2f c and 2f s −3f c . The majority of fault frequencies fell to the average of the spectral lines, and also were not prominent.
The amplitudes of characteristic vibration frequencies in friction torque increase slightly following increased motor speed, based on the theoretical analysis in Section 2.6. However, after modulating to current, the extent of increase is negligible based on the above spectrum analysis.
In fact, from 15 Hz (Figure 16d), the fault frequencies begin to decline. Till 20 Hz (Figure 16e), the fault frequencies are already infrequent and unobvious. The reason is that the ideal motion trail just can be met for lower motor speed. When the motor speed is higher, the sharp collision appears between the ball with raceways. The ideal motion trail of the ball is disturbed by the irregular collision, the feature of multiple characteristic vibration frequencies is also weakened.
Consequently, when the power supply frequency is set below 15 Hz, this is a reasonable option for tracking multiple characteristic vibration frequencies modulated in the current.
In practice, anywhere below 15 Hz is a feasible observation point, provided that there is no sideband overlap with power supply frequency and its harmonics.
Accordingly, the power supply frequency was set to 10 Hz in following investigation.

Spall Size Estimation
Based on the conclusion proposed in Section 2.6, a regular rise and fall should be present in the amplitude of multiple characteristic vibration frequencies in friction torque. Corresponding to the power supply frequency, 10 Hz, the motor speed was focused at 300 r/min. In that case, the theoretical analysis result between spall size and multiple characteristic vibration frequencies in friction torque is as shown in Figure 7b.
After the fluctuant friction torque is modulated to current, the amplitudes of multiple characteristic vibration frequencies will be weakened partly. According to the theoretical analysis result in Section 2.7, the modulation index will decrease gradually, followed with an increase of order of characteristic vibration frequencies. That is to say, the amplitudes of sidebands of 3f c and 2f c will be lower than the result shown in Figure 7b; meanwhile, the amplitude of sidebands of f c will be higher.
In fact, this is verified by the above investigation. In Figure 16b, the bearing with spall length of 10 • , the sideband of 3f c was maximal, but was not very high, like that shown in Figure 7b. Meanwhile, the sideband of f c was not very low, like that shown in Figure 7b. Similar situations are also present in Figure 16a, Figure 16c, and Figure 16d.
Using the experimental platform in Figure 13, the motor was in no-load, and the voltage frequency ratio was always set as 19V/10Hz during the following experiment. Three SKF 6206 bearings were artificially drilled, with an approximate angle of spall length of 5 • , 10 • and 15 • respectively, as shown in Figure 10. A normal bearing, and the above three defective bearings, were installed in the test motor successively. The measured motor speeds were 294 r/min, 293 r/min, 290 r/min, 289 r/min, and the f c values were 17.4 Hz, 17.4 Hz, 17.2Hz, 17.1Hz, accordingly. The squared envelope spectra of the stator currents in all of above cases were as shown as Figure 17.
frequency ratio was always set as 19V/10Hz during the following experiment. Three SKF 6206 bearings were artificially drilled, with an approximate angle of spall length of 5°, 10° and 15° respectively, as shown in Figure 10. A normal bearing, and the above three defective bearings, were installed in the test motor successively. The measured motor speeds were 294 r/min, 293 r/min, 290 r/min, 289 r/min, and the fc values were 17.4 Hz, 17.4 Hz, 17.2Hz, 17.1Hz, accordingly. The squared envelope spectra of the stator currents in all of above cases were as shown as Figure 17.
It should be noted that those frequencies, below power supply frequency, will not be useful as a diagnosis index, because there may be several adjacent frequencies modulated by multiple harmonics present in this area.
In Figure 17a, for a normal bearing, except for fs−fr, there is no fault frequency present in spectrum.
In Figure 17b, for a defective bearing with spall length of 5°, several fault frequencies arise, far beyond the average of spectral lines. The sidebands of fc include: 2fs+fc, 6fs−fc, 4fs+fc and 8fs−fc. There are not obvious sidebands for 2fc and 3fc. One reason is that two impulse responses caused by spall partly overlap each other. Another reason is that the modulation indices of 2fc and 3fc are less than that of fc. In Figure 17c, for defective bearings with spall lengths of 10°, a number of fault frequencies arose. The sidebands of 3fc include: 4fs−3fc and 12fs−3fc. The sidebands of 2fc include: 2fs+2fc and 10fs−2fc. The sidebands of fc include: 5fs−fc, 2fs+fc and 4fs+fc.
According to the trend shown in Figure 7b, as the amplitudes of 3fc and 2fc increased gradually from 5° to 10°, so, the sidebands of 3fc and 2fc emerged increasingly in the current. Because the  It should be noted that those frequencies, below power supply frequency, will not be useful as a diagnosis index, because there may be several adjacent frequencies modulated by multiple harmonics present in this area.
In Figure 17a, for a normal bearing, except for f s −f r , there is no fault frequency present in spectrum. In Figure 17b, for a defective bearing with spall length of 5 • , several fault frequencies arise, far beyond the average of spectral lines. The sidebands of f c include: 2f s +f c , 6f s −f c , 4f s +f c and 8f s −f c . There are not obvious sidebands for 2f c and 3f c . One reason is that two impulse responses caused by spall partly overlap each other. Another reason is that the modulation indices of 2f c and 3f c are less than that of f c .
In Figure 17c, for defective bearings with spall lengths of 10 • , a number of fault frequencies arose. The sidebands of 3f c include: 4f s −3f c and 12f s −3f c . The sidebands of 2f c include: 2f s +2f c and 10f s −2f c .
The sidebands of f c include: 5f s −f c , 2f s +f c and 4f s +f c .
According to the trend shown in Figure 7b, as the amplitudes of 3f c and 2f c increased gradually from 5 • to 10 • , so, the sidebands of 3f c and 2f c emerged increasingly in the current. Because the increase of 3f c is more obvious than others, 4f s −3f c is maximal among all these sidebands; next is 5f s −f c and 2f s +2f c , which is also generally consistent with the trends shown in Figure 7b, considering the additional effect of the modulation index.
In Figure 17d, for defective bearings with spall length of 15 • , there were more fault frequencies present in the spectrum. The sidebands of 2f c include: 2f s +2f c , 6f s −2f c , 2f s −2f c , 10f s −2f c , 4f s +2f c and 8f s −2f c . The sidebands of f c include: 5f s −f c , 2f s +f c and 4f s +f c .
According to the trend shown in Figure 7b, as amplitudes of 3f c are decreased gradually from 10 • to 15 • , so the sidebands of 3f c disappeared in the current. As the amplitudes of 2f c increased slightly from 10 • to 15 • , so, a number of sidebands of 2f c became present in the current. On the other hand, as the amplitudes of f c always increase slightly from 5 • to 15 • based on Figure 7b, so, the 5f s −f c was maximal among all these sidebands; next is 2f s +2f c .
From Figure 17b to Figure 17d, the sidebands of f c are always present, but the amplitude of sidebands did not increase visibly. For this reason, the spall size isn't estimated based on the amplitudes of sidebands of f c .
In general terms, the regular rise and fall of sidebands of 3f c and 2f c in the current are verified based on Figure 17, which is generally consistent with the trend shown in Figure 7b, considering the additional effect of the modulation index. Consequently, some remarkable results can be summarized:

1.
When the spall is a small size (spall length less than 5 • , below a 2.5mm spall length for a 6206 bearing), only the sidebands of f c are visible in the current.

2.
When the spall is a little larger (angle of spall length about 6 • -10 • , approximately a 3-5mm spall length for a 6206 bearing) the sidebands of 2f c and 3f c are gradually more obvious. Following the increase of the spall length, the increase of the sidebands of 3f c becomes more obvious than others; finally, the sidebands of 3f c reach a maximum.

3.
When the spall is large (angle of spall length about 11 • -15 • , approximately a 5.5-7.5mm spall length for a 6206 bearing), following the increase of the spall length, the sidebands of 3f c weaken gradually; finally, the sidebands of 3f c disappear entirely. Meanwhile, the sidebands of 2f c become more obvious gradually.
Consequently, based on the proposed approaches, the spall size in a bearing's outer raceway can be estimated approximately by tracking the multiple characteristic vibration frequencies in the current.
When the spall length enlarges further, the vibration and harsh sound caused by the spall are very obvious; therefore, the larger spall length isn't necessary to investigate.

Conclusions
This paper proposes an approximate fault size estimation approach based on current, which is immune to external environmental interference, addressing the shortcomings of existing approaches which are not suitable for applications in some special and harsh environments. The present research reveals the relation between the spall size and the sidebands modulated in current, proposes a novel idea to estimate spall size by tracking the current sidebands, and introduces a fault-feature-highlighting approach based on reduced voltage frequency ratio. Some noteworthy conclusions are presented as follows: 1.
Compared to traditional fault-feature-highlighting approaches based on instantaneous symmetrical components and Park's vector, the proposed approach has showed better performance.

2.
Below lower power supply frequencies, the proposed fault size estimation approach can effectively track fault size. 3.
The proposed fault size estimation approach only needs to calculate the squared envelope spectrum, an algorithm which is very simple, Therefore, this approach is easy to implement for existing inverter-driven induction motors without complicated calculations and additional sensors, and is suitable for harsh conditions. 4.
Compared to existing fault size estimation approaches, the drawback of the proposed approach is lower estimation accuracy. Therefore, this approach is only suitable for applications in bearing performance degradation assessment and maintenance decision making for induction motors in service.