Acceleration Harmonic Estimation for Hydraulic Servo Shaking Table Based on Multi-Innovation Stochastic Gradient Algorithm

As one of the critical test equipment, hydraulic servo shaking table is widely used in shaking environment simulation of structural components and systems. However, inherent nonlinear factors of a hydraulic servo shaking table can cause amplitude attenuation and phase lag when corresponding to a sinusoidal acceleration signal, which leads to serious harmonic distortion. In order to improve waveform reproduction performance of sinusoidal signals, the amplitude and phase of harmonic should be estimated accurately. In this paper, the multi-innovation stochastic gradient (MISG) algorithm is presented for dynamically estimating the harmonic information. Simulation and experiment results demonstrate that the proposed algorithm has high estimation precision and good convergence performance.


Introduction
Hydraulic servo system plays an indispensable role in industry fields because of its excellent power density, high precision, and fast response, which has attracted the attention from both industry and academic [1][2][3]. As the most important application of the hydraulic servo system, hydraulic servo shaking table is widely used to simulate real-world motion and force within the laboratory. Compared with the field test, laboratorybased structural test is an integral part of the product development process in many fields, such as architecture engineering [4], automobile industry [5], and aircraft [6]. However, the inherent nonlinearities of a hydraulic servo shaking table [7][8][9], such as nonlinear pressure-flow characteristic of servo valve, nonlinear friction of hydraulic actuators, and backlash nonlinearity of mechanical connection, result in serious harmonic distortion when the shaking table is excited by sinusoidal acceleration signals. erefore, it is of great importance to estimate the harmonic information of the hydraulic servo shaking table in a quick and accuracy way.
Various methods were presented to estimate the parameters of harmonic [10], e.g., well-known fast Fourier transform (FFT), least mean square algorithm (LMS), recursive least square (RLS), and Kalman filter algorithm. FFT is widely used in harmonic estimation for its efficiency and simplicity among these methods. However, the nonsynchronous sampling may result in aliasing, spectrum leakage, and picket fence effort [11]. To overcome the above drawbacks, windowed interpolation FFT is employed to improve the harmonic estimation precision, but it cannot provide a satisfactory estimation when suffered from spectrum leakage and harmonic interference [12]. As recursive algorithms, LMS and RLS have several limitations in terms of accuracy and convergence for solving the harmonic estimation problem. e LMS algorithm suffers from poor convergence and being failure in the case of signal drifting and changing conditions [13]. e initialization parameter for the RLS algorithm has a great impact on time varying dynamic signals [14]. e Kalman filter algorithm is capable enough to estimate harmonic parameter in presence of noise and other nonlinearities presented in the response signal [15]. However, the requirement of prior information about the harmonic signal's statistics characteristic and the state matrix's initialization restricts its broader practicable applications [16].
Besides, many hybridized algorithm-based metaheuristic technique and classical optimization method have been reported to estimate the parameters of harmonics [14]. Because the harmonic is linear in amplitude and nonlinear in phase, the hybridized algorithm using the metaheuristic technique estimates the phase, whereas the classical optimization method estimates the amplitude [17]. Liu et al. developed real-time particle swarm optimization (RT-PSO) estimation technique to simultaneously identify all parameters of fundamental signals from measurement data, and then the estimated fundamental signal was used as the reference for harmonics cancellation [18]. Mishra presented a hybrid least square-fuzzy bacterial foraging strategy harmonic estimation algorithm, whereas Takagi-Sugeno fuzzy scheme was adopted to improve the adaptation of the basic foraging strategy [19]. Biswas et al. proposed the combined artificial bee colony-least square algorithm to solve the harmonic estimation problem. Simulation results demonstrated that the proposed algorithm was effective in estimating the magnitude and phase of each harmonic even if the number of harmonics presented in the signal was not known [20]. However, the convergence speed of above proposed algorithms is affected by the initial parameters, and the computational complexity is relatively large, so the estimation accuracy and convergence speed cannot be guaranteed.
Since the gradient estimation method only needs to compute the first-order derivation, the computation load is relatively small. Aiming at the problem of low estimation precision, many improved gradient algorithms have been proposed. Andrei developed an adaptive conjugate gradient algorithm to compute the optimization value of objective function based on the current point [21]. Deng and Wan devised a conjugate gradient algorithm to solve unconstrained optimization problems, and the search direction at each iteration is determined by rectifying the steepest descent direction with the difference between the current iterative points and that between the gradients [22]. Patrinos and Bemporad proposed a dual fast gradientprojection method for solving the embedded linear model predictive control problem subjected to general polyhedral constraints on inputs and states [23]. Although abovemodified algorithms can improve the convergence rate and the estimation accuracy, the computation load is still heavy. For improving the performance of the stochastic gradient algorithm for system identification, the multi-innovation theory has been proposed, whose basic idea is extending the scalar innovation into an innovation vector [24]. Li et al. illustrated the modified MISG algorithm for the identification of Wiener Hammerstein systems with unknown orders linear subsystems and backlash [25]. Ding successively put forward the MISG algorithm, multi-innovation projection algorithm, and interval-varying MISG identification algorithm to solve the parameter estimation problem [26]. e rest of this paper is organized as follows. Section 2 introduces the hydraulic servo shanking table system. Section 3 deduces the MISG algorithm and gives the harmonic estimation scheme based on the MISG algorithm. Simulation results and experiments are depicted in Section 4. Section 5 draws the main points and contributions. Figure 1 shows the experimental setup of the hydraulic servo shake table, which consists of a shaking table and real-time control system. e control system is controlled by a host-target type, where a commercial PC based on the Windows XP system is provided as the host computer, and the target computer is an industrial computer. e control procedure is programmed by Matlab/Simulink and compiled by Microsoft Visual Studio. NET on the host computer and then downloaded to the target computer through the Ethernet adapter using the TCP/IP protocol for real-time execution [27,28]. e signal conditioner compares the control signal and feedback signal and drives the shaking table to generate the desired movement.

e Dynamic Model of Hydraulic Actuator.
e schematic diagram of a hydraulic actuator is shown in Figure 2, which can be simplified as a symmetrical hydraulic cylinder controlled by the servo valve. e dynamic model is based on some fundamental assumptions that the servo valve dynamics are sufficiently fast to be neglected, oil temperature and bulk modulus are constants, and leakage as well as pressure losses of pipeline are negligible. e load flow of the servo valve is a function of spool displacement and load pressure. Linearized load flow equation can be expressed as where q L is the load flow, p L � p 1 − p 2 is the load pressure, and K q and K c are the flow gain and flow pressure coefficient, which can be described as where C d is the coefficient of a valve port, w is the area gradient of valve orifices, ρ is the density of hydraulic oil, and sgn(·) is the sign function. e load flow in a hydraulic cylinder is related to leakage chamber volume variance and compressibility. e flow continuity equation of an actuator can be deduced as where V t � V 1 + V 2 is the total volume, β e is the effective oil bulk modulus, and C t � C i + 0.5C e is the total leakage. Taking inertia force and vicious damping force into consideration, the force balance equation at the piston and specimen is established as follows: Employing Laplace transform to equations (1), (3), and (4), the transfer function from the input x v to output x p can be derived as where K ce � K c + C t is the total flow pressure coefficient and K h � 4β e A 2 /V t is the hydraulic spring stiffness.
Generally speaking, the A 2 /K ce in a characteristic equation is bigger than B, so K ce B/A 2 is far less than 1. en, equation (5) can be simplified as

Mathematical Problems in Engineering
where ω h and ξ h are angular frequency and damping ratio of the electrohydraulic servo system and K v � (K q /A) is the open-loop gain of electrohydraulic servo system. e three-variable controller (TVC) is important to a hydraulic servo shaking table, where the three variables correspond to displacement, velocity, and acceleration. As shown in Figure 3, TVC includes input filter, TVC feedforward controller, and TVC feedback controller. e input filter is utilized to convert the reference acceleration signal to the desired displacement and generate TVC feedforward variable. e TVC feedforward controller is used to extend the frequency bandwidth of the closed loop by adjusting three feedforward parameters K ar , K vr , and K dr . e displacement and acceleration signal are directly obtained by sensors, while the velocity feedback signal is synthesized using a low pass filter with the displacement and a high pass filter with the acceleration. e TVC feedback controller aims at increasing the damping ratio and original frequency hydraulic by adjusting three feedback parameters K af , K vf , and K df .
e TVC was designed based on the position of closed loop of the hydraulic servo system, but the motion of the hydraulic servo shaking table is given by acceleration data. So, the input filter was used to synthesize the desired position drive signal. e transfer function of an input filter is written as follows: where ω 0 � 6.28 and ξ 0 � 0.5 are the initial frequency and damping of acceleration control. According to the zero-poles placement principle, the detail parameters of the TVC controller are determined as follows [29]:

e Sinusoidal Shaking Test.
For a hydraulic servo shaking table, it is quite difficult to replicate the input excitation signal because of the nonlinearities in the system. As a result of nonlinear phenomena, harmonic distortion occurs in the system acceleration response when the hydraulic servo shaking table corresponds to a sine vibration signal. e acceleration response of the hydraulic servo shaking table excited by (4 sin(2π × 5t)m/s 2 ) is shown in Figure 4. e time-domain results demonstrate that the actual response is seriously distorted owing to the nonlinearities. It is easy to see from the frequency domain that the frequencies of harmonics are integer multiples of the fundamental frequency.
at is to say, the second harmonic is located in 10 Hz, the third harmonic in 10 Hz, and so on. As a critical indicator, total harmonic distortion (THD) is used to measure the waveform reproduction performance of sinusoidal signals [30]. e THD analysis result is shown in Table 1. It can be seen that the THD value is 23.4%, which means serious harmonic distortion. e fundamental response is close to the input excitation amplitude. Among the seventh harmonics, the fifth harmonic is the largest harmonic, and the seventh harmonic is the smallest harmonic: where A 1 is the amplitude of fundamental response, A 2 is the amplitude of second harmonic, and A N is the amplitude of nth harmonic.

The Harmonic Estimation Scheme Based on MISG Algorithm
Supposing that the harmonic identification problem can be expressed as the celebrated controlled autoregressive model where y (k) represents the system output, u (k) represents the system input, v (k) represents the stochastic noise with zero mean, d represents the pure delay, and A(z − 1 ) and B(z − 1 ) are polynomials defined in the following form: Define the information vector φ(k) and the parameter vector θ as follows: a 2 , . . . , a n a , b 0 , b 1 , . . . , b n Combining equations (9)- (11), the linear regression model for harmonic estimation of a hydraulic servo shaking table can be deduced as follows: For convenience, supposing that k is the current time, y (k) and φ(k) are the current data, and y(k − i) and φ(k − i) are the past data. e cost function can be expressed as According to the stochastic gradient search principle, the recursive stochastic gradient identification algorithm can be expressed as where θ(k) represents the estimation of θ(k), (1/r(k)) represents the convergence factor, and e(k) represents the innovation. e stochastic gradient algorithm does not involve the calculation of the covariance matrix, which greatly reduces the calculation load. However, the stochastic gradient algorithm only uses the current data, so the convergence speed and precision cannot be guaranteed. In order to overcome the above drawbacks, it is necessary to increase the data utilization on each calculation. From the view of innovation modification, the scalar innovation can be extended into the  innovation vector, which is called the multi-innovation identification algorithm. e proposed identification algorithm uses not only the current data but also the past data [31]. e information vector Φ(p, k), output vector Y(p, k), and innovative vector with the innovation length p can be expressed as In practical applications, the convergence speed and accuracy of the stochastic gradient algorithm are not very ideal. erefore, the forgetting factor (FF) is added to the stochastic gradient algorithm to improve the convergence speed and accuracy. From equations (14)- (16), the FFMISG algorithm with innovation length p can be derived as follows: where λ represents the forgetting factor. For the MISG algorithm, the cost function can be expressed as Compared with the stochastic algorithm only using the current data, the MISG uses not only the current data but also the past data, which improves the convergence performance for the parameter estimation. e MISG algorithm repeatedly utilizes the available data in the neighboring iteration to enhance the estimation accuracy, which is the essential characteristic of the MISG algorithm. e more the length of multi-innovation, the more higher the parameter estimation precision for the same measurement data. In other words, a large multi-innovation length results in a better model accuracy, but the price is a larger computation load [32][33][34][35].

3.2.
e Harmonic Estimation Model. Various nonlinear factors, such as friction, dead zone, and backlash in the hydraulic servo shaking table, can lead to sinusoidal acceleration response distortion when the shaking table is excited by sinusoidal signals, which seriously affects the waveform reproduction performance. In order to suppress harmonic distortion, it is necessary to estimate the amplitude and phase of each harmonic quickly and accurately. In general, the frequency of higher harmonic generated in the sinusoidal vibration test is an integral multiple of the fundamental frequency. us, the sinusoidal acceleration response of the hydraulic servo shaking table can be expressed as where N is the number of harmonics, f is the frequency, and A n (t) and φ n (t) are the amplitude and phase of each harmonic. e discrete time version of sinusoidal acceleration response can be written as e estimated model is constructed as follows: where A n (k) and φ n (k) are the estimated amplitude and phase of A n (k) and φ n (k). Equation (23) can be rewritten as Define , and w n2 (k) � A n (k)sin φ n (k); the information vector x(k) and parameter vector w(k) are given as follows: e difference between the actual value and the estimation value is e(k) � a(k) − a(k). (26) en, the cost function can be expressed as Taking the first-order derivation of J(w) with respect to w gives e FFMISG algorithm-based acceleration harmonic estimation scheme can be illustrated in Figure 5. e reference signal is, respectively, fed to the shaking table and information vector. e desired and estimated acceleration response are compared, and the weight vector is updated through the FFMISG algorithm. e acceleration response will be accurately estimated when the cost function converges to zero.
According to the gradient searching, the MISG parameter estimation algorithm is listed as follows: In this way, amplitude and phase of each harmonic as well as the estimated harmonic are derived as follows:

Simulation Results.
In order to verify the effectiveness of the proposed estimation algorithm, the distorted waveform containing the eighth sine waveform is considered as a test signal. e fundamental frequency is 3 Hz, and the frequency of second to eighth harmonics is integral of the fundamental frequency. e amplitude and phase information of the test signal are given in Table 2.
In the simulation, the sampling period is 1 ms, the forgetting factor is 0.98, and the innovation length is 2. e amplitude and phase estimation results of the test signal using the FFMISG algorithm are shown in Figures 6 and 7, respectively. It can be seen that the estimated values for both amplitude and phase are closer with the actual value in Table 1. In addition, the estimated values of amplitude and phase all converge to the nominal values within 0.5 s.
is shows that the acceleration harmonic estimation scheme is very accurate and efficient. e reconstructed signal-based estimated results and the original test signal are shown in Figure 8(a). Except less deviation in the case of initial stage, both estimated and actual signal matches in the plot. Figure 8(b) shows that the estimation error only exists at the initial stage and quickly converges to zero after 1 s.

Experiment Results.
During the experiment, the proposed estimation algorithm is embedded in the real-time controller in Figure 1.
e estimation algorithm and the basic controller run in the target computer. e acceleration response waveform and the estimation results can be directly achieved from a host computer. In addition, the parameters of FFMISG are kept as same as the simulation. From Figure 4 and Table 3, it can be known that there are seven higher harmonics besides the fundamental when the excitation is (4 sin(2π × 5t)m/s 2 ). Figures 9 and 10 show a comparative estimation of amplitudes and phase of experimental signal using SG and FFMISG algorithms, respectively. From these figures, it is verified that FFMISG estimates more accurately

Mathematical Problems in Engineering 7
Amplitude (   as compared to SG. Moreover, the estimated amplitude results in Figure 9 are consistent well with the values obtained by FFT in Table 3, which further validates the effectiveness of the proposed algorithm. Figure 11 shows the waveform of each harmonic reconstructed from the estimated weight vector. It is easy to see that each harmonic is well estimated after a large fluctuation, which is critical for harmonic suppression. e estimated error is shown in Figure 12. It can be seen from the time domain that residual error stays around zero after a period of time. It can be seen from the frequency-domain response that each harmonic of acceleration response is well counteracted. is means that the estimated signals not only have the right amplitudes but also have the right phases.
us, the proposed estimation algorithm based on FFMISG can dynamically track the harmonics.

Conclusions
Sinusoidal shaking test of equipment using hydraulic shaking tables is a common experimental practice. However, the nonlinear factor of a hydraulic servo shaking table can lead to a higher harmonic appeared in the acceleration response, which causes harmonic response and reduces the waveform reproduction performance. In order to cancel the harmonic distortion, fast and accurate harmonic estimation algorithm is required. To this end, harmonic estimation algorithm based on MISG is developed herein to estimate the amplitude and phase information of higher harmonic.
In order to improve the estimation accuracy and the convergence rate of the SG parameter estimation algorithm, the MISG algorithm is presented by means of expanding the scalar innovation into the innovation vector. In contrast to the SG algorithm, the MISG algorithm can improve the parameter estimation accuracy and promote the convergence during the recursive process. In particular, the proposed harmonic estimation algorithm can be easily combined with the original control scheme and run in the target computer.
is ensures the real-time performance, which is very important for the practical application on a shaking table.
Compared with the FFT result calculated offline, the harmonic information can be calculated online by the proposed MISG algorithm. Since higher harmonics are in less dominance in all harmonics, the harmonic estimation here only considers up to the eighth harmonic including fundamental. Simulation and experimental results revealed that the estimation scheme can not only estimate the amplitude and phase of each harmonic but also reconstruct individual harmonic at the same time.
e obtained identification results can be used as the basis for harmonic cancellation of a hydraulic servo shaking table, which is helpful for the sinusoidal vibration test.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.