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

Authors: Xiancheng Wang; Wei Li; Jianjun Yao; Zhenshuai Wan
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

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.


---

## Body

## 1. 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–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, laboratory-based 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–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. Therefore, 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. The LMS algorithm suffers from poor convergence and being failure in the case of signal drifting and changing conditions [13]. The initialization parameter for the RLS algorithm has a great impact on time varying dynamic signals [14]. The 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 gradient-projection method for solving the embedded linear model predictive control problem subjected to general polyhedral constraints on inputs and states [23]. Although above-modified 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].The rest of this paper is organized as follows. Section2 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.


## 2. The Hydraulic Servo Shaking Table

Figure1 shows the experimental setup of the hydraulic servo shake table, which consists of a shaking table and real-time control system. The 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. The 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]. The signal conditioner compares the control signal and feedback signal and drives the shaking table to generate the desired movement.Figure 1
Experimental setup of the hydraulic servo shaking table.


### 2.1. The Dynamic Model of Hydraulic Actuator

The schematic diagram of a hydraulic actuator is shown in Figure2, which can be simplified as a symmetrical hydraulic cylinder controlled by the servo valve. The 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.Figure 2
Schematic diagram of hydraulic actuator controlled by servo valve.The load flow of the servo valve is a function of spool displacement and load pressure. Linearized load flow equation can be expressed as(1)qL=Kqxv−KcpL,where qL is the load flow, pL=p1−p2 is the load pressure, and Kq and Kc are the flow gain and flow pressure coefficient, which can be described as(2)Kq=∂qL∂xv=Cdwps−sgnxvpLρ,Kc=−∂qL∂pL=Cdwxv14ρps−sgnxvpL,where Cd 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.The load flow in a hydraulic cylinder is related to leakage chamber volume variance and compressibility. The flow continuity equation of an actuator can be deduced as(3)qL=Adxpdtx˙v+Vt4βedpLdt+CtpL,where Vt=V1+V2 is the total volume, βe is the effective oil bulk modulus, and Ct=Ci+0.5Ce 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:(4)ApL=md2xpdt2+Bdxpdt.Employing Laplace transform to equations (1), (3), and (4), the transfer function from the input xv to output xp can be derived as(5)Gs=Kq/AVtm/4βeA2s3+Kcem/A2+VtB/4βeA2s2+KceB/A2+1s,where Kce=Kc+Ct is the total flow pressure coefficient and Kh=4βeA2/Vt is the hydraulic spring stiffness.Generally speaking, theA2/Kce in a characteristic equation is bigger than B, so KceB/A2 is far less than 1. Then, equation (5) can be simplified as(6)Gs=Kq/AsVt/4βeA2s2+mKce/A2+BVt/4βeA2s2+1=Kvss2/ωh2+2ξh/ωhs+1,where ωh and ξh are angular frequency and damping ratio of the electrohydraulic servo system and Kv=Kq/A is the open-loop gain of electrohydraulic servo system.The 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 Figure3, TVC includes input filter, TVC feedforward controller, and TVC feedback controller. The input filter is utilized to convert the reference acceleration signal to the desired displacement and generate TVC feedforward variable. The TVC feedforward controller is used to extend the frequency bandwidth of the closed loop by adjusting three feedforward parameters Kar, Kvr, and Kdr. The 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. The TVC feedback controller aims at increasing the damping ratio and original frequency hydraulic by adjusting three feedback parameters Kaf, Kvf, and Kdf.Figure 3
TVC control schematic of a hydraulic servo shake table.The 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. The transfer function of an input filter is written as follows:(7)Gins=Kus2+v0s+d0=Kud0s2/ω02+2ξ0/ω0s+1,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]:(8)Kdf=ωrωnc2Kvωh2,Kvf=2Kdfξncωnc+Kdfωr−1Kv,Kaf=2Kdfξncωrωnc+Kdfωnc2−2ξhKvωhKdr=2Kdfξncωnc,Kvr=Kdf,Kar=2Kdfωnc2.


### 2.2. The 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.The acceleration response of the hydraulic servo shaking table excited by4sin2π×5tm/s2 is shown in Figure 4. The 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. That 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]. The THD analysis result is shown in Table 1. It can be seen that the THD value is 23.4%, which means serious harmonic distortion. The 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:(9)THD=A12+A22+⋯+AN2A1×100%,where A1 is the amplitude of fundamental response, A2 is the amplitude of second harmonic, and AN is the amplitude of nth harmonic.Figure 4
Acceleration response of the sinusoidal shaking test: (a) time-domain response; (b) frequency spectrum response.
(a)
(b)Table 1
THD analysis results of the sinusoidal shaking test.
THD
Harmonica amplitude (m/s2)
23.4%
Fundamental
Second
Third
Fourth
Fifth
Sixth
Seventh
Eighth
3.873
0.365
0.279
0.222
0.640
0.367
0.085
0.102


## 2.1. The Dynamic Model of Hydraulic Actuator

The schematic diagram of a hydraulic actuator is shown in Figure2, which can be simplified as a symmetrical hydraulic cylinder controlled by the servo valve. The 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.Figure 2
Schematic diagram of hydraulic actuator controlled by servo valve.The load flow of the servo valve is a function of spool displacement and load pressure. Linearized load flow equation can be expressed as(1)qL=Kqxv−KcpL,where qL is the load flow, pL=p1−p2 is the load pressure, and Kq and Kc are the flow gain and flow pressure coefficient, which can be described as(2)Kq=∂qL∂xv=Cdwps−sgnxvpLρ,Kc=−∂qL∂pL=Cdwxv14ρps−sgnxvpL,where Cd 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.The load flow in a hydraulic cylinder is related to leakage chamber volume variance and compressibility. The flow continuity equation of an actuator can be deduced as(3)qL=Adxpdtx˙v+Vt4βedpLdt+CtpL,where Vt=V1+V2 is the total volume, βe is the effective oil bulk modulus, and Ct=Ci+0.5Ce 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:(4)ApL=md2xpdt2+Bdxpdt.Employing Laplace transform to equations (1), (3), and (4), the transfer function from the input xv to output xp can be derived as(5)Gs=Kq/AVtm/4βeA2s3+Kcem/A2+VtB/4βeA2s2+KceB/A2+1s,where Kce=Kc+Ct is the total flow pressure coefficient and Kh=4βeA2/Vt is the hydraulic spring stiffness.Generally speaking, theA2/Kce in a characteristic equation is bigger than B, so KceB/A2 is far less than 1. Then, equation (5) can be simplified as(6)Gs=Kq/AsVt/4βeA2s2+mKce/A2+BVt/4βeA2s2+1=Kvss2/ωh2+2ξh/ωhs+1,where ωh and ξh are angular frequency and damping ratio of the electrohydraulic servo system and Kv=Kq/A is the open-loop gain of electrohydraulic servo system.The 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 Figure3, TVC includes input filter, TVC feedforward controller, and TVC feedback controller. The input filter is utilized to convert the reference acceleration signal to the desired displacement and generate TVC feedforward variable. The TVC feedforward controller is used to extend the frequency bandwidth of the closed loop by adjusting three feedforward parameters Kar, Kvr, and Kdr. The 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. The TVC feedback controller aims at increasing the damping ratio and original frequency hydraulic by adjusting three feedback parameters Kaf, Kvf, and Kdf.Figure 3
TVC control schematic of a hydraulic servo shake table.The 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. The transfer function of an input filter is written as follows:(7)Gins=Kus2+v0s+d0=Kud0s2/ω02+2ξ0/ω0s+1,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]:(8)Kdf=ωrωnc2Kvωh2,Kvf=2Kdfξncωnc+Kdfωr−1Kv,Kaf=2Kdfξncωrωnc+Kdfωnc2−2ξhKvωhKdr=2Kdfξncωnc,Kvr=Kdf,Kar=2Kdfωnc2.


## 2.2. The 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.The acceleration response of the hydraulic servo shaking table excited by4sin2π×5tm/s2 is shown in Figure 4. The 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. That 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]. The THD analysis result is shown in Table 1. It can be seen that the THD value is 23.4%, which means serious harmonic distortion. The 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:(9)THD=A12+A22+⋯+AN2A1×100%,where A1 is the amplitude of fundamental response, A2 is the amplitude of second harmonic, and AN is the amplitude of nth harmonic.Figure 4
Acceleration response of the sinusoidal shaking test: (a) time-domain response; (b) frequency spectrum response.
(a)
(b)Table 1
THD analysis results of the sinusoidal shaking test.
THD
Harmonica amplitude (m/s2)
23.4%
Fundamental
Second
Third
Fourth
Fifth
Sixth
Seventh
Eighth
3.873
0.365
0.279
0.222
0.640
0.367
0.085
0.102


## 3. The Harmonic Estimation Scheme Based on MISG Algorithm

### 3.1. The MISG Algorithm

Supposing that the harmonic identification problem can be expressed as the celebrated controlled autoregressive model(10)Az−1yk=z−dBz−1uk+vk,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 Az−1 and Bz−1 are polynomials defined in the following form:(11)Az−1na=1+a1z−1+⋯+anz−na,Bz−1nb=b0++b1z−1+⋯+bnz−nbnb0≠0,nb≤na.Define the information vectorφk and the parameter vector θ as follows:(12)φk=−yk−1,…,−yk−na,uk−d,…,uk−d−nbT,θ=a1,a2,…,ana,b0,b1,…,bnb.Combining equations (9)–(11), the linear regression model for harmonic estimation of a hydraulic servo shaking table can be deduced as follows:(13)yk=φTkθ+vk.For convenience, supposing thatk is the current time, y (k) and φk are the current data, and yk−i and φk−i are the past data. The cost function can be expressed as(14)Jθ=Eyk−φkθ2.According to the stochastic gradient search principle, the recursive stochastic gradient identification algorithm can be expressed as(15)θ^k=θ^k−1+φkrkek,(16)ek=yk−φTkθ^k−1,(17)rk=rk−1+φk2,where θ^k represents the estimation of θk, 1/rk represents the convergence factor, and ek represents the innovation.The 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. The proposed identification algorithm uses not only the current data but also the past data [31].The information vectorΦp,k, output vector Yp,k, and innovative vector with the innovation length p can be expressed as(18)Φp,k=φk,φk−1,…,φk−p+1,Yp,k=yk,yk−1,…,yk−p+1T,Ep,k=Yp,k−ΦTp,kθ^k−1.In practical applications, the convergence speed and accuracy of the stochastic gradient algorithm are not very ideal. Therefore, 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:(19)θ^k=θ^k−1+Φp,krkEp,k,Ep,k=Yp,k−ΦTp,kθ^k−1,rk=λrk−1+Φp,k2,where λ represents the forgetting factor.For the MISG algorithm, the cost function can be expressed as(20)Jθ=EYp,k−ΦTp,kθ2.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. The 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. The 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–35].


### 3.2. The 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. Thus, the sinusoidal acceleration response of the hydraulic servo shaking table can be expressed as(21)at=∑n=1NAntsin2πnft+φnt,where N is the number of harmonics, f is the frequency, and Ant and φnt are the amplitude and phase of each harmonic.The discrete time version of sinusoidal acceleration response can be written as(22)ak=∑n=1NAnksin2πnfk+φnk.The estimated model is constructed as follows:(23)a^k=∑n=1NA^nksin2πnfk+φ^nk,where A^nk and φ^nk are the estimated amplitude and phase of Ank and φnk. Equation (23) can be rewritten as(24)ak=∑n=1NAnkcosφnksin2πnfk+Anksinφnkcos2πnfk.Definexn1k=sin2πfk, xn2k=cos2πfk, wn1k=Ankcosφnk, and wn2k=Anksinφnk; the information vector xk and parameter vector wk are given as follows:(25)xk=x11k,x12k,…,xN1k,xN2kT,wk=w11k,w12k,…,wN1k,wN2kT.The difference between the actual value and the estimation value is(26)ek=ak−a^k.Then, the cost function can be expressed as(27)Jw=12e2k.Taking the first-order derivation ofJw with respect to w gives(28)gradJw=∂Jw∂w=∂Jw∂w11,∂Jw∂w12,…,∂Jw∂wN1,∂Jw∂wN2.The FFMISG algorithm-based acceleration harmonic estimation scheme can be illustrated in Figure5. The reference signal is, respectively, fed to the shaking table and information vector. The desired and estimated acceleration response are compared, and the weight vector is updated through the FFMISG algorithm. The acceleration response will be accurately estimated when the cost function converges to zero.Figure 5
Structure of harmonic estimation based on the FFMISG algorithm.According to the gradient searching, the MISG parameter estimation algorithm is listed as follows:(29)w^k=w^k−1+Φp,krkEp,k,Ep,k=ek,ek−1,…,ek−p+1T,ek−i=yk−i−∑j=1nxjTkw^jk,i=0,1,…,p−1,Φp,k=xk,xk−1,…,xk−p+1,rk=λrk−1+Φp,k2.In this way, amplitude and phase of each harmonic as well as the estimated harmonic are derived as follows:(30)A^nk=w^n12k+w^n22k,φ^nk=tan−1w^n2kw^n1k,a^nk=xn1kw^n1k+xn2w^n2k.


## 3.1. The MISG Algorithm

Supposing that the harmonic identification problem can be expressed as the celebrated controlled autoregressive model(10)Az−1yk=z−dBz−1uk+vk,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 Az−1 and Bz−1 are polynomials defined in the following form:(11)Az−1na=1+a1z−1+⋯+anz−na,Bz−1nb=b0++b1z−1+⋯+bnz−nbnb0≠0,nb≤na.Define the information vectorφk and the parameter vector θ as follows:(12)φk=−yk−1,…,−yk−na,uk−d,…,uk−d−nbT,θ=a1,a2,…,ana,b0,b1,…,bnb.Combining equations (9)–(11), the linear regression model for harmonic estimation of a hydraulic servo shaking table can be deduced as follows:(13)yk=φTkθ+vk.For convenience, supposing thatk is the current time, y (k) and φk are the current data, and yk−i and φk−i are the past data. The cost function can be expressed as(14)Jθ=Eyk−φkθ2.According to the stochastic gradient search principle, the recursive stochastic gradient identification algorithm can be expressed as(15)θ^k=θ^k−1+φkrkek,(16)ek=yk−φTkθ^k−1,(17)rk=rk−1+φk2,where θ^k represents the estimation of θk, 1/rk represents the convergence factor, and ek represents the innovation.The 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. The proposed identification algorithm uses not only the current data but also the past data [31].The information vectorΦp,k, output vector Yp,k, and innovative vector with the innovation length p can be expressed as(18)Φp,k=φk,φk−1,…,φk−p+1,Yp,k=yk,yk−1,…,yk−p+1T,Ep,k=Yp,k−ΦTp,kθ^k−1.In practical applications, the convergence speed and accuracy of the stochastic gradient algorithm are not very ideal. Therefore, 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:(19)θ^k=θ^k−1+Φp,krkEp,k,Ep,k=Yp,k−ΦTp,kθ^k−1,rk=λrk−1+Φp,k2,where λ represents the forgetting factor.For the MISG algorithm, the cost function can be expressed as(20)Jθ=EYp,k−ΦTp,kθ2.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. The 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. The 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–35].


## 3.2. The 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. Thus, the sinusoidal acceleration response of the hydraulic servo shaking table can be expressed as(21)at=∑n=1NAntsin2πnft+φnt,where N is the number of harmonics, f is the frequency, and Ant and φnt are the amplitude and phase of each harmonic.The discrete time version of sinusoidal acceleration response can be written as(22)ak=∑n=1NAnksin2πnfk+φnk.The estimated model is constructed as follows:(23)a^k=∑n=1NA^nksin2πnfk+φ^nk,where A^nk and φ^nk are the estimated amplitude and phase of Ank and φnk. Equation (23) can be rewritten as(24)ak=∑n=1NAnkcosφnksin2πnfk+Anksinφnkcos2πnfk.Definexn1k=sin2πfk, xn2k=cos2πfk, wn1k=Ankcosφnk, and wn2k=Anksinφnk; the information vector xk and parameter vector wk are given as follows:(25)xk=x11k,x12k,…,xN1k,xN2kT,wk=w11k,w12k,…,wN1k,wN2kT.The difference between the actual value and the estimation value is(26)ek=ak−a^k.Then, the cost function can be expressed as(27)Jw=12e2k.Taking the first-order derivation ofJw with respect to w gives(28)gradJw=∂Jw∂w=∂Jw∂w11,∂Jw∂w12,…,∂Jw∂wN1,∂Jw∂wN2.The FFMISG algorithm-based acceleration harmonic estimation scheme can be illustrated in Figure5. The reference signal is, respectively, fed to the shaking table and information vector. The desired and estimated acceleration response are compared, and the weight vector is updated through the FFMISG algorithm. The acceleration response will be accurately estimated when the cost function converges to zero.Figure 5
Structure of harmonic estimation based on the FFMISG algorithm.According to the gradient searching, the MISG parameter estimation algorithm is listed as follows:(29)w^k=w^k−1+Φp,krkEp,k,Ep,k=ek,ek−1,…,ek−p+1T,ek−i=yk−i−∑j=1nxjTkw^jk,i=0,1,…,p−1,Φp,k=xk,xk−1,…,xk−p+1,rk=λrk−1+Φp,k2.In this way, amplitude and phase of each harmonic as well as the estimated harmonic are derived as follows:(30)A^nk=w^n12k+w^n22k,φ^nk=tan−1w^n2kw^n1k,a^nk=xn1kw^n1k+xn2w^n2k.


## 4. Simulation and Experiment Results

### 4.1. 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. The fundamental frequency is 3 Hz, and the frequency of second to eighth harmonics is integral of the fundamental frequency. The amplitude and phase information of the test signal are given in Table2.Table 2
Harmonic information of the test signal.
Harmonic order
Frequency (Hz)
Amplitude (m/s2)
Phase (rad)
Fundamental
3
8
−0.1
Second harmonic
6
7
−0.7
Third harmonic
9
6
0.5
Fourth harmonic
12
5
−0.6
Fifth harmonic
15
4
1.2
Sixth harmonic
18
3
−0.4
Seventh harmonic
21
2
0.7
Eighth harmonic
24
1
−0.8In the simulation, the sampling period is 1 ms, the forgetting factor is 0.98, and the innovation length is 2. The amplitude and phase estimation results of the test signal using the FFMISG algorithm are shown in Figures6 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. This shows that the acceleration harmonic estimation scheme is very accurate and efficient.Figure 6
Estimated amplitude of the test signal.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure 7
Estimated phase of the test signal.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)The reconstructed signal-based estimated results and the original test signal are shown in Figure8(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.Figure 8
Simulation results: (a) estimated signal; (b) estimated error.
(a)
(b)


### 4.2. Experiment Results

During the experiment, the proposed estimation algorithm is embedded in the real-time controller in Figure1. The estimation algorithm and the basic controller run in the target computer. The 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 4sin2π×5tm/s2. 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 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.Table 3
Main parameters of a hydraulic servo shaking table.
Parameter
Value
K
q
1.45 × 10−3 (m3·s−1)/V
K
c
2 × 10−12 (m3/s·Pa)
A
1.88 × 10−3 m2
B
2.5 × 104 N/(m·s−1)
β
e
6.9 × 108 Pa
V
t
9.6 × 10−4 m3
m
100 kg
K
af
1.6 × 10−4
K
vf
8 × 10−6
K
df
6 × 10−2
K
ar
5.6 × 10−5
K
vr
6 × 10−8
K
ar
7.2 × 10−2
K
u
1.2Figure 9
Estimated amplitude of harmonics.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure 10
Estimated phases of harmonics.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure11 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. The 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. This means that the estimated signals not only have the right amplitudes but also have the right phases. Thus, the proposed estimation algorithm based on FFMISG can dynamically track the harmonics.Figure 11
Estimated waveform of harmonics.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure 12
Estimation error of the sinusoidal shaking test: (a) time-domain error response; (b) frequency spectrum response error.
(a)
(b)


## 4.1. 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. The fundamental frequency is 3 Hz, and the frequency of second to eighth harmonics is integral of the fundamental frequency. The amplitude and phase information of the test signal are given in Table2.Table 2
Harmonic information of the test signal.
Harmonic order
Frequency (Hz)
Amplitude (m/s2)
Phase (rad)
Fundamental
3
8
−0.1
Second harmonic
6
7
−0.7
Third harmonic
9
6
0.5
Fourth harmonic
12
5
−0.6
Fifth harmonic
15
4
1.2
Sixth harmonic
18
3
−0.4
Seventh harmonic
21
2
0.7
Eighth harmonic
24
1
−0.8In the simulation, the sampling period is 1 ms, the forgetting factor is 0.98, and the innovation length is 2. The amplitude and phase estimation results of the test signal using the FFMISG algorithm are shown in Figures6 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. This shows that the acceleration harmonic estimation scheme is very accurate and efficient.Figure 6
Estimated amplitude of the test signal.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure 7
Estimated phase of the test signal.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)The reconstructed signal-based estimated results and the original test signal are shown in Figure8(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.Figure 8
Simulation results: (a) estimated signal; (b) estimated error.
(a)
(b)


## 4.2. Experiment Results

During the experiment, the proposed estimation algorithm is embedded in the real-time controller in Figure1. The estimation algorithm and the basic controller run in the target computer. The 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 4sin2π×5tm/s2. 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 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.Table 3
Main parameters of a hydraulic servo shaking table.
Parameter
Value
K
q
1.45 × 10−3 (m3·s−1)/V
K
c
2 × 10−12 (m3/s·Pa)
A
1.88 × 10−3 m2
B
2.5 × 104 N/(m·s−1)
β
e
6.9 × 108 Pa
V
t
9.6 × 10−4 m3
m
100 kg
K
af
1.6 × 10−4
K
vf
8 × 10−6
K
df
6 × 10−2
K
ar
5.6 × 10−5
K
vr
6 × 10−8
K
ar
7.2 × 10−2
K
u
1.2Figure 9
Estimated amplitude of harmonics.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure 10
Estimated phases of harmonics.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure11 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. The 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. This means that the estimated signals not only have the right amplitudes but also have the right phases. Thus, the proposed estimation algorithm based on FFMISG can dynamically track the harmonics.Figure 11
Estimated waveform of harmonics.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure 12
Estimation error of the sinusoidal shaking test: (a) time-domain error response; (b) frequency spectrum response error.
(a)
(b)


## 5. 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. This 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. The 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.



---

*Source: 3063469-2020-03-17.xml*