Crack Monitoring in Rotating Shaft Using Rotational Speed Sensor-Based Torsional Stiffness Estimation with Adaptive Extended Kalman Filters

Abstract

In this study, we present an alternative solution for detecting crack damages in rotating shafts under torque fluctuation by directly estimating the reduction in torsional shaft stiffness using the adaptive extended Kalman filter (AEKF) algorithm. A dynamic system model of a rotating shaft for designing AEKF was derived and implemented. An AEKF with a forgetting factor (λ) update was then designed to effectively estimate the time-varying parameter (torsional shaft stiffness) owing to cracks. Both simulation and experimental results demonstrated that the proposed estimation method could not only estimate the decrease in stiffness caused by a crack, but also quantitatively evaluate the fatigue crack growth by directly estimating the shaft torsional stiffness. Another advantage of the proposed approach is that it uses only two cost-effective rotational speed sensors and can be readily implemented in structural health monitoring systems of rotating machinery.

1. Introduction

Rotating machinery (or turbomachinery) has steadily been in the field of interest for industrial applications in internal combustion engines, power generators, turbines, and high-speed machining [1]. Rotating machinery generally consists of a rotor and a non-rotating part (stator), with torque transmitted through a rotating shaft. Cracks in rotary shafts are among the most dangerous and significant defects. The crack occurs in rotating shafts because of various mechanisms such as high and low cycle fatigue, stress corrosion, or unbalanced force caused by the rotor offset [2]. The shafts of the above-mentioned machines are typically subjected to harsh working conditions, such as loading and temperature variations. Thus, successive failures can lead to enormous economic and human resource losses. If a crack propagates continuously and is not detected in advance, an abrupt failure may occur, leading to catastrophic consequences. Thus, real-time monitoring of crack damage in the rotating shaft is essential.

Generally, contact sensors, which provide high data accuracy and convenience, can be used to detect such cracks in a rotating shaft. However, rotating and internal parts are generally difficult to measure directly. Thus, it is challenging to monitor shaft cracks using a contact sensor, such as a strain gauge. Therefore, in recent years, fault diagnosis studies on rotating machinery have focused on indirect detection methods through vibration response characteristic analysis of components, such as bearings and gears. As a result, numerous vibration-based crack detection techniques have been developed over the last decades [3]. These techniques include experimental signal-based and model-based methods. Several model-based crack detection methods, such as wavelet transform [4], have been developed to enhance fault diagnosis. Experimental signal-based methods using nonlinear vibration responses have also been widely used for damage detection in structures [5,6].

For more precise, reliable, and effective detection, various non-destructive techniques (NDT), such as radiography, magnetic particle inspection, and ultrasonic methods, are attempted to diagnose and monitor the behavior of rotating machines, although these techniques consume more time and are expensive [7]. However, high-frequency amplitudes are too small to detect cracks, and responses can be generated by assembly tolerances, manufacturing state noise, and other defects.

These disadvantages of the current technology necessitate developing non-traditional technology for detecting structural surface damage, such as cracks in the rotating shaft [8,9,10,11,12]. Non-model-based crack detection has also been attempted as a statistical-based data analysis method, using trained models from artificial neural networks [13] and genetic algorithms [14]. Recently, machine learning has been studied as a solution for detecting defects effectively without human experts [15]. However, this method is data-inefficient because we cannot acquire sufficient experimental data on actual crack sequences for large systems. Over the last decades, some research on structural health monitoring has been conducted based on the adaptive extended Kalman filter algorithm (AEKF) [16,17,18,19,20]. For example, the Kalman filter with the forgetting factor method had been applied to several systems, such as a lithium-ion battery, to consider the variation of system model parameters [21]. However, it is still necessary to study a new detection method, although previous studies have shown promising results in detecting cracks in rotating shafts.

Therefore, this study primarily aims to provide an alternative solution for detecting crack damages in rotating shafts by directly estimating the change in stiffness using the adaptive extended Kalman filter algorithm (AEKF) with a forgetting factor update. To the best of our knowledge, we report for the first time that it is possible to achieve a new means of detecting the torsional crack in a rotating shaft using AEKF with a forgetting factor update algorithm. Cracks of varying geometry are caused by different types of stress-field directions and are classified according to their orientation with respect to the shaft axis, as shown in Figure 1a. The direction of the stress field depends on the type of stress (such as bending or torsion) and geometric factors. When high cyclic stress is repeated, the crack propagates such that the crack plane is perpendicular to the direction of the tensile stress field. When bending stress is applied to the shaft, a stress field forms along the axis, and the crack propagates into the shaft section, creating a transverse crack, which is frequently called a breathing crack [22]. Torsional stress forms a tensile stress field in the direction of 45° to the shaft axis. In this study, we focused on torsional slant cracks of shafts and attempted to use Kalman-filter-based torsional stiffness estimation. When fatigue cracks occur in rotating shaft systems under alternating torque excitation, the cracks gradually grow larger over time as they are repeatedly opened and closed. As the cross-sectional area decreases, the torsional stiffness of the shaft suddenly decreases, as shown in Figure 1b. The AEKF-based estimator of shaft torsional stiffness using a dynamic model of rotating machinery is described in Section 2. Simulation results using the proposed algorithm under sinusoidal torque input are presented in Section 3. The simulation results were experimentally validated, as described in Section 4.

2. Design of Adaptive Extended Kalman Filters

2.1. Dynamic Modeling of Rotating Shaft

Although there are various connecting structures, such as bearings and shafts, in a real rotating shaft system, the system model was formulated based on a dynamic circular shaft, to which torque and rotational speed were applied together. It is assumed to be a lumped-parameter model with two primary masses (i.e., a semi-definite system with a single natural frequency) because it is unlikely to excite the higher modes in our simple test bed system (single frequency excitation). The system elements simulated the driving-load motor dynamo for experimental verification of the proposed algorithm. The stiffness of the bellows coupling connecting the shaft and the damping effect of the bearing stand were neglected. As shown in Figure 2, the rotating shaft model comprised four components. The governing equation for the shaft rotation is given as follows:

$$J_m{\ddot{\theta }}_m+c_m{\dot{\theta }}_m+k_s({\theta }_m-{\theta }_l)=T_m,$$

(1)

$$k_s({\theta }_m-{\theta }_l)-J_l{\ddot{\theta }}_l=0,$$

(2)

where $J_m$ is the moment of inertia (driving motor), $J_l$ is the moment of inertia (load motor), $k_s$ is the torsional stiffness of the shaft, and $c_m$ is the damping coefficient of the viscous friction of the driving motor. When torque is applied by the load motor, the angular velocity difference between the two sides is caused by the stiffness of the shaft connecting the two motors. It is necessary to express the dynamic model into a state-space model to implement the Kalman filter algorithm.

$$\dot{x}=Ax+Bu,$$

(3)

$$y=Cx+Du,$$

(4)

where x is the state vector, u is the input vector, and $\dot{x}$ is the time derivative of the state vector. In Equation (3), A is a state matrix, and B is an input matrix. In Equation (4), y is a measurement variable, C is a measurement matrix, and D is a feed-forward matrix.

The four state variables for stiffness estimation are selected as shown in Equation (5). In this study, the time-varying stiffness is treated as a state variable, and it is assumed to be linearly proportional to the crack sizes for designing Kalman filters, although it can be changed by the nonlinear dynamics of the rotating shaft system [23]. ${\theta }_m-{\theta }_l$ is the difference in angular displacement on both sides and ${\omega }_l$ is the angular velocity of the load motor. The fourth state variable is the angular velocity of the driving motor. The driving motor torque is an input for the system. From Equations (1) and (2), the state space equation was derived as follows:

$$x={\left[x_1{x}_2{x}_3{x}_4\right]}^T={\left[\left({\theta }_m-{\theta }_l\right){\omega }_l{k}_s{\omega }_m\right]}^T,$$

(5)

$$u=T_m,$$

(6)

$$\dot{x}=f(x,u)=\left[\begin{array}{c}{\dot{x}}_1\\ {\dot{x}}_2\\ {\dot{x}}_3\\ {\dot{x}}_4\end{array}\right]=\left[\begin{array}{c}{x}_4-x_2\\ \frac{x_3{x}_1}{J_l}\\ 0\\ \frac{-c_m{x}_4-x_3{x}_1+u}{J_m}\end{array}\right],$$

(7)

$$y=h(x)=\left[\begin{array}{c}{x}_2\\ x_4\end{array}\right]=\left[\begin{array}{cccc}\begin{array}{l}0\\ 0\end{array}& \begin{array}{l}1\\ 0\end{array}& \begin{array}{l}0\\ 0\end{array}& \begin{array}{l}0\\ 1\end{array}\end{array}\right]\left[\begin{array}{c}{\theta }_m-{\theta }_l\\ {\omega }_l\\ k_s\\ {\omega }_m\end{array}\right].$$

(8)

The reformulated system model $f(x,u)$ is a nonlinear model, and the measurement model h(x) is a linear model with an actual measurable angular velocity value as the output. For the estimation of torsional stiffness, an extended Kalman filter (EKF) that linearizes a nonlinear model is required, and an adaptive EKF (AEKF) with a P-adaptive loop is proposed to improve estimation performance. As the proposed AEKF algorithm is based on the discrete-time domain, the continuous equation was discretized using the Euler method, as shown in Equation (9).

$$\dot{x}=\frac{x(k)-x(k-1)}{\Delta t} \to x(k)=\dot{x}(k)\Delta t+x(k-1),$$

(9)

where Δt is the time step, and k and k − 1 represent the time instant at t = kΔt and t = (k − 1)Δt, respectively. Substituting Equation (7) into Equation (9), Equations (10) and (11) are defined as follows:

$$\left\{\begin{array}{l}{x}_k=f_{k-1}(x_{k-1},u_{k-1})\\ y_k=h_k(x_k)\end{array}\right..$$

(10)

2.2. Adaptive Extended Kalman Filters

Kalman filtering is a state-estimation technique developed by Rudolf Kalman in 1960. It features a recursive structure and optimally estimates the state of a linear dynamic system based on measurements contaminated by noises. Kalman filters are used in many industrial fields, such as computer vision, robotics, and vehicular electronics [24,25]. The general linear discrete-time system model required to design the KF is given as

$$\begin{array}{l}{x}_{k+1}=A{x}_k+B{u}_k+w_k\\ y_k=H{x}_k+v_k\end{array},$$

(11)

where $w_k$ is a multivariate Gaussian distribution system noise variable with a covariance matrix, and $v_k$ is a multivariate Gaussian distribution measurement noise variable with a covariance matrix. In this study, there was no input in the measurement model, and the application of the EKF was based on the nonlinear model. The general discrete-time equation is as follows:

$$\left\{\begin{array}{l}{x}_k=f_{k-1}(x_{k-1},u_{k-1})+w_{k-1}\\ y_k=h_k(x_k)+v_k\end{array}\right..$$

(12)

The extended Kalman filter assumes differentiability of the state-change function instead of linearity of the model. The nonlinear system model was linearized using the Jacobian, and the Jacobian matrix was calculated based on the previous estimate.

$$A_{k-1}={\left.\frac{\partial f_{k-1}}{\partial x}\right|}_{{\widehat{x}}_{k-1}}{B}_{k-1}={\left.\frac{\partial f_{k-1}}{\partial u}\right|}_{{\widehat{x}}_{k-1}},$$

(13)

$$H(k)={\left.\frac{\partial h_k}{\partial x}\right|}_{{\widehat{x}}_{k|k-1}}.$$

(14)

Matrices $A$ and $B$ of the rotating shaft system model linearized using Equations (13) and (14) are as follows:

$$A^{*}=\frac{\partial f(x,u)}{\partial x}=\left[\begin{array}{cccc}0& -1& 0& 1\\ \frac{x_3}{J_L}& 0& \frac{x_1}{J_L}& 0\\ 0& 0& 0& 0\\ -\frac{x_3}{J_m}& 0& -\frac{x_1}{J_m}& \frac{C_m}{J_m}\end{array}\right],$$

(15)

$$B^{*}=\frac{\partial f(x,u)}{\partial u}={\left[\begin{array}{cccc}0& 0& 0& \frac{1}{J_m}\end{array}\right]}^T,$$

(16)

$$H=\left[\begin{array}{l}\begin{array}{cccc}0& 1& 0& 0\end{array}\\ \begin{array}{cccc}0& 0& 0& 1\end{array}\end{array}\right],$$

(17)

where (*) indicates the system model matrix linearized using the Jacobian.

The discrete-time EKF algorithm has the following form:

▪

Initial estimation stage at $k=0$

$$\left\{\begin{array}{l}{\hat{x}}_0=E[x_0]\\ P_0=E[(x_0-{\hat{x}}_0){(x_0-{\widehat{x}}_0)}^T]\end{array}\right..$$

(18)

where E represents the expected value of the random variable.

▪

Prediction stage

$$\begin{array}{l}\widehat{x}(k|k-1)=f_{k-1}({\widehat{x}}_{k-1},u_{k-1},0)\\ P(k|k-1)=A(k,k-1)P(k-1)A^T(k,k-1)+Q(k-1)\end{array}.$$

(19)

A matrix was differentiated using the Jacobian in Equation (13). In the prediction stage, the variables predicted are a priori state variable and an error covariance matrix.

▪

Correction stage

$$K(k)=P(k|k-1)H^T(k){(H(k)P(k|k-1)H^T(k)+R(k))}^{-1},$$

(20)

$$\begin{array}{l}\widehat{x}(k)=\widehat{x}(k|k-1)+K(k)[y(k)-h_k({\widehat{x}}_{k|k-1},u_k,0)]\\ P(k)=[I-K(k)H(k)]P(k|k-1)\end{array}.$$

(21)

In general, it is difficult to estimate the time-varying parameter (shaft stiffness) using the EKF because filter estimation relies on past data, and state estimation can diverge when past data are not adequate for recursive estimation methods. In this study, an AEKF with a forgetting factor (λ) was used to resolve this technical limitation [26]. The updated forgetting factor corrects the error covariance matrix, and the Kalman gain matrix is increased by the inverse of the forgetting factor. In general, the forgetting factor is considered a constant tuning parameter. However, convergence decreases when the uncertainty is large, such as in a nonlinear model. In this study, an adaptive loop was employed for more weighting to recent data using the residual between the measured and estimated values [27,28]. The AEKF equation is identical to the EKF in Equation (19), except for the forgetting factor in the error covariance equation.

$$P(k+1|k)=\lambda (k+1)A(k+1,k)P(k)A^T(k+1,k)+Q(k),$$

(22)

with $\lambda (k)\ge 1$. Thus, divergence is prevented by considering the influence of the most recently measured data on the state and parameter. The performance of the AEKF is the most important factor because it completely depends on the forgetting factor. The residual $z(k)$ is defined as the difference between the measured and predicted values of the measurement. The residual is a white noise sequence when the optimal filtering gain is used.

$$z(k)=y(k)-H(k)\widehat{x}(k|k-1).$$

(23)

For any gain, the covariance of the residuals is expressed as:

$$C_0(k)=E[z(k)z^T(k)]=H(k)P(k|k-1)H^T(k)+R(k).$$

(24)

The auto covariance of the residual is

$$\begin{array}{l}{C}_j(k)=E[z(k+j)z^T(k)]\\ \hspace{1em} =H(k+j)A(k+j,k+j-1)\\ \hspace{1em}\hspace{1em}\hspace{1em} \times [I-K(k+j-1)H(k+j-1)]\cdots A(k+2,k+1)\\ \hspace{1em}\hspace{1em}\hspace{1em} \times [I-K(k+1)H(k+1)]A(k+1,k)\\ \hspace{1em}\hspace{1em}\hspace{1em} \times [P(k|k-1)H^T(k)-K(k)C_0(k)]\\ \forall j=1,2,3,\cdots \end{array}.$$

(25)

In general, $C_j(k)$ in Equation (25) is equal to zero when Equations (20) and (24) are substituted into Equation (25), implying that the residual sequences are uncorrelated when the optimal gain is applied. However, the actual covariance of the residual $C_0(k)$ is different from the theoretical covariance, owing to errors in the system model parameters and noise covariance. Therefore, $C_j(k)$ may not be equal to zero. In Equation (25), we can choose a forgetting factor such that the last term of $C_j(k)$ for all is zero.

$$P(k|k-1)H^T(k)-K(k)C_0(k)=0.$$

(26)

In the optimal condition, $S(k)$ and $g(\lambda ,k)$ are as follows:

$$S(k)=P(k|k-1)H^T(k)-K(k)C_0(k),$$

(27)

$$g(\lambda ,k)=\frac{1}{2}{\displaystyle \sum _{i=1}^n{\displaystyle \sum _{j=1}^m{S}_{ij}^2(k)}}.$$

(28)

The optimality of the Kalman filter can be determined through Equation (27), which is a scalar function, and $S_{ij}(k)$ is $(i,j)$ th element of $S(k)$. As the smaller $g(k)$ yields more optimal filter, the forgetting factor $\lambda (k)$ should be selected to minimize $g(k)$.

Various studies have been conducted based on the least-squares estimation (LSE) approach to better track time-varying parameters of dynamic systems. In this study, a recursive estimation method with a forgetting factor update was introduced to track time-varying parameters. The constant forgetting factor was optimally updated based on the following equation (i.e., the gradient descent method):

$${\lambda }^{l+1}(k)={\lambda }^l(k)+\phi \frac{\partial g^l(\lambda ,k)}{\partial {\lambda }^l(k)}\forall l=0,1,2,\dots ,$$

(29)

with initial conditions

$${\lambda }^0(1)=1,{\lambda }^0(k)=\lambda (k-1),$$

(30)

where k is the time series and l is the iteration time of the time instant. φ is the step length (i.e., learning rate, $0<\phi <1$). If Equation (31) is satisfied in the p-th iteration (i.e., converges), the iteration is stopped, and the optimal forgetting factor is determined using Equation (32).

$$\left|{\lambda }^p(k)-{\lambda }^{p-1}(k)\right|<\epsilon .$$

(31)

$$\lambda (k)=\max \{1,{\lambda }^p(k)\}.$$

(32)

However, the iterative numerical method does not guarantee real-time processing. Finally, a one-step AEFK algorithm was used to resolve this computational burden. In the given system, state Equations (12), (18), and (19) have the following assumption:

Assumption 1. 

$Q(k),R(k)$ and $P(0)$ are positive definite.

Assumption 2. 

The measurement matrix$H(k)$is fully ranked, and the optimal forgetting factor can be calculated as

$$\lambda (k)=\max \left\{1,trace[N(k)]/trace[M(k)]\right\},$$

(33)

where

$$M(k)=H(k)A(k,k-1)P(k-1)A^T(k,k-1)H^T(k),$$

(34)

$$N(k)=C_0(k)-H(k)Q(k-1)H^T(k)-R(k).$$

(35)

The $C_0(k)$ value was estimated using the recursive equation as follows:

$$C_0(k)=G_1(k)/G_2(k),$$

(36)

$$G_1(k)=G_1(k-1)/\lambda (k-1)+z(k)z^T(k),$$

(37)

$$G_2(k)=G_2(k-1)/\lambda (k-1)+1$$

(38)

with initial conditions $G_1(0)=0$ and $G_2(0)=0$. The proofs of Equations (33)–(35) was derived by substituting Equation (20), which derives the Kalman gain value into Equation (26).

$$\begin{array}{l}P(k|k-1)H^T(k)\\ \times \{I-{[H(k)P(k|k-1)H^T(k)+R(k)]}^{-1}{C}_0(k)\}=0\end{array}$$

(39)

$$H(k)P(k|k-1)H^T(k)=C_0(k)-R(k).$$

(40)

Equation (40) implies that, with Assumptions 1 and 2, the optimality condition described in Equation (26) is equivalent to Equation (24). Substituting Equation (14) into Equation (40), and then reconstructing it yields the following:

$$\begin{array}{l}\lambda (k)H(k)A(k,k-1)P(k-1)A^T(k,k-1)H^T(k)\\ =C_0(k)-H(k)Q(k-1)H^T(k)-R(k)\end{array}.$$

(41)

The overall estimation process using the AEKF algorithm with a forgetting factor update is shown in Figure 3.

3. Estimation of Shaft Torsional Stiffness

3.1. Simulation Scheme

Based on the proposed algorithm, a situation in which cracks occur owing to shaft damage was simulated using MATLAB^®. The parameter values of the system model are listed in

Table 1. Parameters for the estimation of torsional stiffness.

Parameters (Unit)	Value
Inertia moment of load motor $J_l(Nm{m}^2)$	580
Inertia moment of driving motor $J_m(Nm{m}^2)$	180
Damping constant cm (Nmm·s/rad)	1000
Shaft torsional stiffness ks (Nmm/rad)	735,000

. In this study, to mimic crack fatigue due to persistent cyclic excitation, a sinusoidal torque input was applied (frequency of 1 Hz, $T_m(t)=\mathrm{10,000}\sin (2\pi t) Nmm$). The angular velocity measurement data from the simulation model was set to be contaminated by the white Gaussian random noise $v(k)=N(0,{({10}^{-3})}^2)$.

To evaluate the response time of the proposed estimator, step response to a sudden downward step input is used for the crack initiation scenario, which corresponds to large cracks in an experiment. Then, it was assumed that the torsional stiffness suddenly decreased from 735,000 to 345,000 Nmm/rad in 10 s, as shown in Figure 4.

The initial state value and error covariance for estimation are as follows:

$$\begin{array}{l}{x}_0=\begin{array}{cccc}[0& 0& \mathrm{735,000}& 0]\end{array}\\ P_0=diag\begin{array}{cccc}([0.01& 1& \mathrm{800,000}& 1])\end{array} \end{array}.$$

(42)

The system noise covariance matrix Q and the measurement noise covariance R for Equations (21) and (22) were tuned in various cases as follows, and the optimal estimates were derived:

$$Q=\left[\begin{array}{cccc}{10}^{-8}& & & \\ & {10}^{-7}& & \\ & & {10}^{-7}& \\ & & & {10}^{-7}\end{array}\right] ,R=\left[\begin{array}{cc}{10}^{-3}& \\ & {10}^{-3}\end{array}\right] .$$

(43)

To evaluate the basic estimation performance of the AEKF, the root-mean-squared error (RMSE) at the kth time instant was calculated for a more rigorous analysis.

$$RMSE(k)=\sqrt{\frac{1}{k}{\displaystyle \sum _{i=1}^k{(p(i)-\widehat{p}(i))}^2}},$$

(44)

where k is the time instant at t = kΔt, $p(i)$ and $\widehat{p}(i)$ are the true (i.e., Figure 4) and estimated values, respectively. The steady-state mean of the RMSE (MRMSE) was then calculated to exclude the effect of transient behavior. The basic estimation results for the sudden torsional stiffness drop are shown in Figure 5. The AEKF accurately estimated the sudden torsional stiffness change. In contrast, the EKF did not track the time-varying shaft stiffness change. The forgetting factor was appropriately changed by the P-adaptive loop when the stiffness rapidly decreased in 10 s. Additional scenarios with different reduction rates are applied to the simulation model to investigate the effectiveness of the proposed algorithm. These different scenarios allow for the evaluation of the tracking performance of the proposed algorithm under the same conditions, such as process and measurement noise covariance matrices. As shown in Figure 6, a gradual reduction from $1.5\times {10}^4 \mathrm{Nmm}/\mathrm{rad}$ starts at approximately 5 s, drops to $0.4\times {10}^4 \mathrm{Nmm}/\mathrm{rad}$ at 35 s (simulating a situation where the crack is propagating). When a crack growth is propagating and a gradual torsional stiffness drop occurs, the AEKF can deal appropriately.

3.2. Robustness Analysis

The robustness of the proposed estimation model under noise and parametric model uncertainty was analyzed by introducing perturbations to sensor noise and main parameters. To evaluate the robustness under noise and parametric uncertainties, the relative error to the nominal value (i.e., normalized performance measure) was quantitatively calculated.

$$\mathrm{Relative} \mathrm{Error}=\frac{\left|MRMS{E}_{\mathrm{Perturbed}}-MRMS{E}_{\mathrm{Nominal}}\right|}{MRMS{E}_{\mathrm{Nominal}}}.$$

(45)

As sensor information is inherently contaminated by electrical noise, the effect of electrical noise on the estimated performance was examined. The sensor data were contaminated by adding a white Gaussian random noise. A probability density function is shown in Figure 7a,b as an example. Considering the random noise (error) distribution can be fitted to a normal Gaussian distribution with variance (${\sigma }^2=0.00197$, Case 1; ${\sigma }^2=0.00298$, Case 2), it was confirmed by white Gaussian random noise. The proposed AEKF appeared to be robust against the Gaussian random noise extracted from the sensor data because the estimation results appeared to be similar to the original data, with no significant discrepancy, as shown in Figure 7c,d.

The estimation performance of the proposed AEKF model was evaluated under parametric uncertainty, such as the moment of inertia. The moment of inertia on both sides is an important model uncertainty because it depends on the size, weight, and connection structure of the coupling. The nominal value for the moment of inertia of the load motor (580 $\mathrm{N}\mathrm{m}{\mathrm{m}}^2$) was perturbed by −20% (464 $\mathrm{N}\mathrm{m}{\mathrm{m}}^2$) and +20% (696 $\mathrm{N}\mathrm{m}{\mathrm{m}}^2$), and the nominal inertia moment of the driving motor (180 $\mathrm{N}\mathrm{m}{\mathrm{m}}^2$) was also perturbed by −20% (144 $\mathrm{N}\mathrm{m}{\mathrm{m}}^2$) and +20% (216 $\mathrm{N}\mathrm{m}{\mathrm{m}}^2$). The damping coefficient varied under normal operating conditions (800~1200 Nmm·s/rad), depending on the bearing lubrication condition. As the estimation results were similar to the nominal values within a reasonable range under various parametric uncertainties, the robustness of the proposed model was demonstrated, as shown in Figure 8.

4. Experimental Validation

4.1. Experimental Set-Up

The proposed shaft health-monitoring method was experimentally validated using a torque dynamo. An aluminum hollow-rod specimen ($D_o$: 20 mm, $D_i$: 18.2 mm, L: 550 mm) was used for the rotating shaft, as shown in Figure 9. The torque dynamo comprises a driving motor and torque-controlled load motor (Mitsubishi HG-SR152, 10 Hz bandwidth). Sinusoidal torque ($T_m(t)=\mathrm{10,000}\sin (2\pi t)$ Nmm) was applied at a rotating speed of 5.23 rad/s (50 RPM). For the shaft crack scenario, the shaft was exchanged in turn from a normal shaft without cracks to a cracked shaft in the 45° direction (Figure 9b). The crack depth was set to 5 mm to ensure that the shaft stiffness could suddenly drop from the original value. To examine the possibility of applying a non-contact angular velocity sensor (tachometer), the measurement model considered the angular velocity values on both sides of the rotating shaft. The angular velocities of both sides were measured using a photoelectric detector-type rotational velocity sensor (ONO SOKKI, model: LG-930), which calculates the rotation speed by counting the light reflected on the gear per rotation as a pulse. The real-time monitoring performance was evaluated using a dSPACE^® system (DS1104).

In this study, the recursive least square estimator (RLSE) was used to identify the unknown model parameters. The rotating shaft model expressed in Equations (1) and (2) were reformulated in the matrix form as follows:

$$y_k=h_k^T{\theta }_k+v_k$$

(46)

where

$$y_k=\left[\begin{array}{c}{T}_m\\ 0\end{array}\right], h_k^T=\left[\begin{array}{llll}{\ddot{\theta }}_m& 0& {\dot{\theta }}_m& \left({\theta }_m-{\theta }_l\right)\\ 0& -{\ddot{\theta }}_l& 0& \left({\theta }_m-{\theta }_l\right)\end{array}\right], {\theta }_k={[J_m, J_l, c_m, k_s]}^T$$

(47)

In addition to the measured data from the sensors, other information was required for the two matrices $y_k$ and $h_k^T$. First, the input torque ($T_m$) in matrix $y_k$ was measured using an in-line torque sensor (model: YDR-2K), as shown in Figure 9. The angular displacement (${\theta }_m-{\theta }_l$) and two angular accelerations (${\ddot{\theta }}_m$, ${\ddot{\theta }}_l$) for the matrix $h_k^T$ was obtained by directly differentiating and integrating using the low-pass filtering of the angular velocity signal. The RLSE was then designed as follows:

▪

Initial estimates

$${\widehat{\theta }}_0=E\left[\theta \right]$$

(48)

$$P_0=E\left[\left(\theta -{\widehat{\theta }}_0\right){\left(\theta -{\widehat{\theta }}_0\right)}^T\right]$$

(49)

▪

Kalman gain calculation

$$K_{k+1}=P_k{h}_{k+1}{\left(h_{k+1}^T{P}_k{h}_{k+1}+w_{k+1}^{-1}\right)}^{-1}$$

(50)

▪

Parameter update

$${\widehat{\theta }}_{k+1}={\widehat{\theta }}_k+K_{k+1}\left(y_{k+1}-h_{k+1}^T{\widehat{\theta }}_k\right)$$

(51)

▪

Covariance update

$$P_{k+1}=\left(I-K_{k+1}{\mathit{h}}_{k+1}^T\right)P_k$$

(52)

All parameters of the rotating system are successfully estimated because they converge to a steady-state final positive value after 5000 iterations, as shown in Figure 10. The identified system parameters of the rotating shaft model are listed in

Table 2. Identified system parameters of the rotating shaft model.

Parameters (Unit)	Value
Inertia moment of load motor $J_l({\mathrm{Nmm}}^2)$	595
Inertia moment of driving motor $J_m({\mathrm{Nmm}}^2)$	20
Damping constant cm (Nmm·s/rad)	280
Shaft torsional stiffness ks (Nmm/rad)	15,000

.

4.2. Results and Discussion

For the AEKF estimation model, the initial states were set, and the two noise covariance matrices (Q and R) were tuned by trial and error, as listed in

Table 3. Tuning parameters for the AEKF estimation model.

P0	diag[0.1 1 650,000 1]
Q	diag[1 2.1 2 2.1] × 10−5
R	diag[9 9] × 10−8

. The shaft stiffness estimated using the proposed algorithm was compared in Figure 11. The estimated stiffness became steady-state and converged after 15 s in both cases. In the case of the normal state (no crack), the convergence value was identical to the system identification value (i.e., 15,000 Nmm/rad). When the stiffness changes owing to the sudden drop of crack (crack depth 5 mm) from 15,000 Nmm/rad (normal) to a certain value (abnormal crack, in this case approximately 7500), the proposed algorithm can detect this sudden drop. However, it was difficult to conform to the shift in shaft stiffness by naked eyes from the two angular velocity inputs.

A different crack scenario was established to further investigate the effectiveness of the proposed algorithm. The crack depth was further increased (11 mm) to suddenly drop from 15,000 Nmm/rad to below 7500 Nmm/rad owing to the reduction in cross-sectional area (the ratio of the crack segment area to the original cross-sectional area was 65%) [29,30]. Similar to Figure 11, the proposed algorithm can track this stiffness drop due to the heavy crack, as shown in Figure 12. The proposed estimation model could not only estimate the decrease in stiffness caused by a crack, but also quantitatively evaluate the fatigue crack growth by directly estimating the shaft torsional stiffness. The robustness of the proposed estimation model under noise uncertainty was evaluated by introducing the perturbations in sensor noise. The original sensor signal was filtered by a digital moving average filter (no phase delay). Two corrupted signals were generated; the low-pass filtering is small (Case 1, less contaminated) and off (Case 2, more contaminated, i.e., raw data). The probability density distribution of sensor noise extracted from the original sensor signal was similar to Gaussian distribution, as shown in Figure 13c. The proposed estimation model seemed to be robust against the Gaussian random noise in all sensor data because the estimation results appeared to be similar regardless of the degree of contamination, as shown in Figure 13a,b. In the case of heavy crack (crack depth 11 mm), the proposed estimation model turned out to be more robust against the Gaussian random noise, as shown in Figure 13d. The robustness of the proposed estimation model under model uncertainty was not investigated in the experiment because it was unlikely to significantly change two main model parameters (inertia moment of load and driving motor).

5. Conclusions

In this study, the torsional crack in the rotating shaft was successfully detected in real-time by estimating the reduction of torsional stiffness in the rotating shaft using the AEKF approach with forgetting factor update. The main contributions of this study are summarized as follows:

▪

We concluded that the proposed approach is a promising alternative means for detecting torsional cracks in rotating shafts despite the difficulty in tuning the Q and R matrices of the AEKF.

▪

The proposed estimation model could not only estimate the decrease in stiffness caused by a crack but also quantitatively evaluate the fatigue crack growth by directly estimating the shaft torsional stiffness.

▪

Another advantage of the proposed approach is that it uses only two cost-effective rotational speed sensors; therefore, it does not require noncontact-type torque sensors, which are typically expensive and suffer from durability limitations.

With these advantages, the proposed approach can be readily implemented in structural health monitoring systems of rotating machinery. In future research, we will continue to address some of the ongoing issues. In particular, the localization of cracks in rotating shafts should be studied further. In addition, if the input variables cannot be measured, an advanced algorithm should be applied to simultaneously estimate unknown input and state variables.