# Nonsmooth Recursive Identification of Sandwich Systems with Backlash-Like Hysteresis

Authors: Ruili Dong; Yonghong Tan; Hui Chen; Yangqiu Xie
Journal: Journal of Applied Mathematics
 (2012)

---

## Abstract

A recursive gradient identification algorithm based on the bundle method for sandwich systems with backlash-like hysteresis is presented in this paper. In this method, a dynamic parameter estimation scheme based on a subgradient is developed to handle the nonsmooth problem caused by the backlash embedded in the system. The search direction of the algorithm is estimated based on the so-called bundle method. Then, the convergence of the algorithm is discussed. After that, simulation results on a nonsmooth sandwich system are presented to validate the proposed estimation algorithm. Finally, the application of the proposed method to anX-Y moving positioning stage is illustrated.


---

## Body

## 1. Introduction

Usually, a sandwich system with backlash-like hysteresis is defined as the system that a backlash-like hysteresis is sandwiched between two linear dynamic subsystems. In engineering applications, many mechanical systems such as mechanical transmission systems, servo control systems, and hydraulic valve systems can be described by the so-called sandwich systems with backlash-like hysteresis. The reason to cause the backlash-like hysteresis phenomenon is mainly due to the gaps existing in transmission mechanism systems such as gearbox and ball screw.Recently, identification of sandwich systems has become one of the interesting issues in the domain of modeling and control for complex systems. References [1–3] proposed the recursive identification methods for the sandwich system with smooth nonlinearities. The main ideas of those approaches are to extend the linear system identification methods to smooth nonlinear cases. Moreover, there have been some methods for the identification of Hammerstein or Wiener systems with backlash-like hysteresis [4–8], most of which are the modified linear system identification methods.However, until today, there have been very few publications concerning the identification of the sandwich systems with backlash-like hysteresis. Reference [9] proposed a method to identify the sandwich systems with backlash-like hysteresis, but the approach is still based on idea to extend the linear system identification method to nonlinear cases. On the other hand, the switching functions in that method have significant influence on the convergence speed of the algorithm.In this paper, a recursive gradient algorithm based on the bundle method is proposed to identify parameters of the sandwich model. In this algorithm, the effect of the nonsmoothness caused by the backlash-like hysteresis in sandwich system is considered. In order to obtain the optimizing search direction at the nonsmooth points of the system, the Clarke subgradient technique is utilized based on the idea of the bundle method [10–12]. By comparing with the above-mentioned available methods, the proposed method employs the nonsmooth optimization technique to identify the nonsmooth sandwich systems with backlash-like hysteresis. Thus, it will provide us with a new approach for dealing with on-line modeling of nonsmooth dynamic systems. A numerical example will be presented to evaluate the performance of the proposed approach. Finally, experimental results on an X-Y moving positioning stage are illustrated.


## 2. Brief Description of Sandwich Systems with Backlash

The structure of a sandwich system with backlash-like hysteresis is shown in Figure1, in which a backlash-like hysteresis is embedded between the input and output linear subsystems, that is, L1(·) and L2(·). It is assumed that input u(k) and output y(k) can be measured directly, but the internal variables x(k) and v(k) are not measurable.Figure 1
The structure of the sandwich system with backlash-like hysteresis.Suppose that both linear subsystems are stable, and the time delaysq1 and q2 in L1(·) and L2(·) are known, respectively. The corresponding discrete-time models of L1(·) and L2(·) are, respectively, written as
(2.1)x(k)=-∑i2=1naai2x(k-i2)+∑j2=0nbbj2u(k-q1-j2),y(k)=-∑i1=1ncci1y(k-i1)+∑j1=0nddj1v(k-q2-j1),
where na and nb are the orders of L1(·), q1 is the time delay, and ai2 as well as bj2 are the coefficients of L1(·); nc and nd are the orders of L2(·), q2 is the time delay, and ci1 and dj1 are the coefficients of L2(·). Let both b0 and d0 be equal to unity for unique representation.Note that the backlash-like hysteresis shown in Figure1 is specified by the slopes m1 and m2 as well as the absolute thresholds, D1 and D2, where 0<m1<∞,0<m2<∞,0<D1<∞, and 0<D2<∞. Hence, the discrete-time model of the backlash-like hysteresis is described as
(2.2)v(k)={m1(x(k)-D1),x(k)>v(k-1)m1+D1,x(k)>x(k-1),increasezone,v(k-1),v(k-1)m2-D2≤x(k)≤v(k-1)m1+D1,memoryzone,m2(x(k)+D2),x(k)<v(k-1)m2-D2,x(k)<x(k-1),decreasezone.For the convenience to describe the system, the discrete-time model of the backlash-like hysteresis can be rewritten as(2.3)m(k)=m1+(m2-m1)g(k),v1(k)=m(k)(x(k)+g(k)x(k)-D1g1(k)+D2g2(k)),v(k)=v1(k)+[v(k-1)-v1(k)](g1(k)-1)(g2(k)-1),
where the switching functionsg(k), g1(k), and g2(k) are, respectively, defined as
(2.4)g(k)={0,Δx(k)>01,Δx(k)≤0,g1(k)={1,x(k)>v(k-1)m1+D1,x(k)>x(k-1),0,else,g2(k)={1,x(k)<v(k-1)m2-D2,x(k)<x(k-1),0,else,
where Δx(k)=x(k)-x(k-1).Thus, (2.1)–(2.3) present the model to describe the sandwich system with backlash-like hysteresis. Hence, the unknown parameter vector of the model can be written as θ∈Rna+nb+nc+nd+4, where
(2.5)θ=[c1,…,cnc,a1,…,anam1,m2,D1,D2,b1,…,bnb,d1,…,dnd]T.
According to concept of the gradient algorithm, define the objective function as
(2.6)Q(k,θ^(k))=∑k=1n[y(k)-y^(k,θ^(k))]22=12∑k=1nf(k,θ^(k)),
where θ^ is the estimate of θ, and y^(k,θ^(k)) is the output of system model. The optimal estimate of θ^ can be obtained by minimizing the above-mentioned criterion.


## 3.  The Nonsmooth Estimation of the Sandwich Model with Backlash-Like Hysteresis

In this section, a gradient-based identification algorithm is proposed for identification of the sandwich system with backlash-like hysteresis. Due to the nonsmoothness of the backlash, the gradients of the system output with respect to the parameters of the backlash at nonsmooth points will not exist. The smooth gradient-based methods directly applied to nonsmooth systems may fail in convergence [13]. On the other hand, the genetic algorithms [14] or Powell’s method [15], which are based on derivative-free techniques, may be unreliable and become inefficient when the system structure is complicated. Thus, we should find a special way for solving this problem. The simplest way to solve the problem is to apply the Clarke subgradients [11] to the approximation of the gradients at the nonsmooth points.The basic idea of the bundle method is to approximate the subdifferential ofQ(k,θ^(k)) with respect to θ^(k) by gathering the subgradients from previous iterations into a bundle for the nonsmooth objective function Q(k,θ^(k)). The gradient ∇Q(k,θ^(k)) can change discontinuously, and some change of the gradient may not be small in the neighborhood of the minimum of the function. So the values of Q(k,θ^(k)) and ∂Q(k,θ^(k)) at a single point θ^(k) do not offer sufficient information of the local behavior of Q(k,θ^(k)). The detail of the bundle method can be found in [10–12] and reference therein.Considering that the sandwich system with backlash-like hysteresis is locally Lipschitz continuous, we have the following definition.Definition 3.1 (see [11]).
LetF: Rn×R→R be locally Lipschitz continuous. This allows one to define a Clarke subgradient of F at ξ as dF(ξ):
(3.1)dF(ξ)∈∂F(ξ),subjectto∂F(ξ)=conv{∇F(ξi)∣ξi⟶ξ,∇F(ξi)exists},
where “conv” denotes the convex hull of a set.
The set of all the Clarke subgradients is the Clarke subdifferential ofFat ξ which is denoted by ∂F(ξ) [11].
Considering that backlash-like hysteresis is a nonsmooth mapping, the gradients of parameters inL1(·) with respect to v(k) do not exist at a nonsmooth point. Hence, we define the parameters of the backlash-like hysteresis and L1(·) as σ={m1,m2,D1,D2,a1⋯ana,b1⋯bnb}∈R4+na+nb. Considering the cost function described by (2.6), the gradients of f(·) with respect to σ will not exist at the nonsmooth points. Hence, at the nonsmooth points of Q(·), the Clarke subdifferential of f(·) with respect to σ, that is, ∂f(σ), can be obtained by
(3.2)∂f(σ)=-conv{[y(k)-y^(k,θ^(k))](∑j1=0nbdj1∂v^(k-j1-q2)(σ))},
where ∂v^(k-j1-q2)(σ)=conv{∇v^(k-j1-q2)(σ)}, and ∇v^(k-j1-q2)(σ) is the gradient of v^, the output of backlash-like hysteresis, with respect to σ at the smooth points. Thus, the corresponding gradients of v^(k-j1-q2) with respect to σ at the smooth points are
(3.3)∇v^(k-j1-q2)(σ)={[u(k-j1-q2-q1)-D^1(k-1),0,-m^1(k-1),0,-m^1(k-1)x^(k-1-j1-q2),…,-m^1(k-1)x^(k-na-j1-q2),m^1(k-1)u(k-1-j1-q2-q1),…,m^1(k-1)u(k-nb-j1-q2-q1)]T,inincreasezones,[0,0,0,0,0,…,0,0,…,0]T,inmemoryzones,[0,u(k-j1-q2-q1)+D^2(k-1),0,m^2(k-1),-m^2(k-1)x^(k-1-j1-q2),…,-m^2(k-1)x^(k-na-j1-q2),m^2(k-1)u(k-1-j1-q2-q1),…,m^2(k-1)u(k-nb-j1-q2-q1)]T,indecreasezones,
where x^(k)=-∑i2=1naa^i2x^(k-i2)+∑j2=0nbb^j2u(k-q1-j2), and the coefficients a^i2 and b^j2 are the corresponding estimated values at the previous step.
Hence, based on (3.2) and (3.3), the Clarke subdifferential of f(·) with respect to σ can be obtained at nonsmooth points of the system. Besides, as L2(·) is a smooth function, the gradients of f(·) with respect to the parameters of the linear subsystems L2(·) always exist. So, the Clarke subdifferential of f(·) with respect to all the unknown parameters of the sandwich system can be determined.
The proper Clarke subgradient directiont(k,θ^(k)) of f(·) with respect to the parameters to be estimated at nonsmooth points can be derived based on
(3.4)min⁡φ,d⁡(φ(k)+12‖t(k,θ^(k))‖2)s.t.-βj(k)+〈hj(k),t(k,θ^(k))〉≤φ(k),∀j∈Jk,
where ∥·∥ denotes the Euclidean norm; Jk is a nonempty subset of {1,…,k}; set φ(k) is the predicted amount of descent; hj(k)∈∂f(k,θ^j*(k)) for j∈Jk, and θ^j*(k) are some trail points (from the past iterations); βj(k)=max{|αj(k)|,γ(sj(k))2} is the locality measure of subgradient; γ≥0 is the distance measure parameter (γ=0 if f(k,θ^(k)) is convex), αj(k)=f(k,θ^(k))-f(k,θ^j*(k))-hj(k)(θ^(k)-θ^j*(k)) is the linearization error; sj(k)=∥θ^j(k)-θ^j*(k)∥+∑i=jk-1∥θ^i+1(k)-θ^i(k)∥ is the distance measure to estimate ∥θ^(k)-θ^j*(k)∥ without the requirement to store the trial point θ^*(k).
According to formula (3.4), t(k) and φ(k) are obtained, that is,
(3.5)t(k,θ^(k))=-∑j∈Jkλjkhj(k)=∑j∈Jk[y(k)-y^(k,θ^j(k))]λjkwj(k)=e(k,θ^(k))h-j(k),(3.6)φ(k)=-‖t(k,θ^(k))‖2-∑j∈Jkλjkβj(k),
where wj(k)=∂y^(k,θ^)/∂θ|θ^=θ^j(k), e(k,θ^(k))=y(k)-y^(k,θ^j(k)), h-j(k)=∑j∈Jkλjkwj(k),λjk≥0, and ∑j∈Jkλjk=1.Remark 3.2.
Iff(·) is convex, the model f(k,θ^(k)) is an underestimate for f(·), and the nonnegative linearization error αj(k) measures the performance of an approximation of the model to the original cost function. If f(·) is nonconvex, these facts are not valid anymore because αj(k) may have a small or even negative value, although the trial point θ^j*(k) locates far away from the current iteration point θ^(k), and thus, the corresponding subgradient hj(k) is worthless. For these reasons, the locality measure of subgradient βj(k) is introduced.
Therefore, the proposed recursive gradient estimation algorithm based on bundle method for the sandwich model with backlash-like hysteresis is shown as follows.Step 1.
Select starting pointθ0∈Rna+nb+nc+nd+4 and stopping parameter δ>0.Calculate f(k,θ0) and vector hj(k)∈∂f(k,θ0), where j∈Jk, Jk={k0}, |Jk|≤k1, |Jk| is the element number of Jk, and k1 is a given positive number. Set βj(k)=0, k=k0 and the line search parameters
(3.7)q∈(0,0.5),q*∈(q,1),η(0)∈(0,1].Step 2.
Calculate optimal solution(φ(k),t(k,θ^(k))) based on formulas (3.2)–(3.6). If φ(k)≥-δ, then stop.Step 3.
Search for the largest step sizeη(k)∈[0,1] such that η(k)≥η(0) and if
(3.8)f(k,θ^(k)+η(k)t(k))≤f(k,θ^(k))+qη(k)φ(k),
it holds
(3.9)φ(k)=f(k,θ^(k)+t(k))-f(k,θ^(k))<0.
Then, we take a long step and set θ^(k+1)=θ^(k)+η(k)t(k) and θ^*(k+1)=θ^(k+1); go to Step 4.
Otherwise, if0<η(k)<η(0), and formula (3.8) holds, then we take a short step and set θ^(k+1)=θ^(k)+η(k)t(k), and θ^*(k+1)=θ^(k)+η*(k)t(k) where η*(k)>η(k). Go to Step 5.
Ifη(k)=0, and formula (3.8) holds, we take a null step, and namely set θ^(k+1)=θ^(k) and θ^*(k+1)=θ^(k)+η*(k)t(k); go to Step 5.Step 4.
LetJk=Jk∪{k+1}, k=k+1; if k≤k1, then Jk={1,…,k}, and if k>k1, then Jk=Jk-1∪{k}∖{k-k1}, then go to Step 2.Step 5.
J
k
=
J
k
∪
{
k
+
1
}, k=k+1; if k≤k1, then Jk={1,⋯,k}, and if k>k1, then Jk=Jk-1∪{k}∖{k-k1}, and the proper Clarke subgradient hj(k) satisfies
(3.10)-βj(k,θ^(k))+hjT(k,θ^(k))t(k,θ^(k-1))≥q*φ(k,θ^(k-1)),
then go to Step 2.Remark 3.3.
In long step, there is an obvious decrease in the value of the objective function. Hence, it is unnecessary to detect discontinuities in the gradient off(·). Thus, we just set hj(k)∈∂f(k,θ^(k)). On the other hand, in short steps and null steps, the gradient of f(·) is discontinuous. Then, based on (3.10), both θ^(k) and θ^*(k) located on the opposite sides of this discontinuity are guaranteed, and the new subgradient hj(k)∈∂f(k,θ*(k)) will force an obvious modification of the next search direction. Hence, the algorithm approximates the effectively searching direction at nonsmooth points based on the bundle method, which cannot be realized by the smooth optimization techniques.Remark 3.4.
If the value ofη(0) is too small, the convergence speed will be very sluggish, while η(0) is too large, and the algorithm may not be convergent. Hence, it is important for η(0) to be chosen properly. Usually, η(0) is chosen based on an empirical method.Remark 3.5.
If all the Clarke subgradients are included inJk, the corresponding storage capacity is infinite. Hence, the number of the subgradients in Jk must be constrained. In the proposed algorithm, we give the upper bound of |Jk|≤k1, and the upper bound k1 is specified by empirical method.


## 4. Convergence of the Estimation

For the convergence of the above-mentioned estimation algorithm, we have the followingTheorem 4.1.
Suppose thatη(k) and βj(k) satisfy
(4.1)0≤η(k)≤2e2(k)h-(k)h-T(k)-βj(k)e2(k)h-T(k)h-(k)(1+[h-(k)h-T(k)]),(4.2)βj(k)≤2e2(k)[h-(k)h-T(k)],
respectively, then the parameters θ can be convergent to a local optimal value.Proof.
The proof of this theorem can be found in Appendix.


## 5. Simulation

The proposed approach is used to identify a numerical sandwich system with backlash-like hysteresis based on the measured system input and output. Suppose that the parameters of the backlash-like hysteresis in the system arem1=1,m2=1.2,D1=0.5, and D2=0.6. The linear subsystems L1(·) and L2(·) are
(5.1)x(k)=-0.1x(k-1)-0.2x(k-2)+1.5u(k-1),(5.2)y(k)=-1.2y(k-1)-0.32y(k-2)+2v(k-1)-0.1v(k-2),
respectively.That impliesa1=0.1,a2=0.2,b0=1.5,c1=1.2,c2=0.32,d0=2, and d1=-0.1. In the simulation, both b0 and d0 are assumed to be equal to unity for model uniqueness, which implies that the corresponding equivalent true values of the coefficients are a~1=0.1, a~2=0.2, c~1=1.2, c~2=0.32, d~1=-0.05, m~1=3, m~2=3.6, D~1=0.33, and D~2=0.4, respectively, but this does not affect the properties of the whole system.In the simulation, the signal to excite the system is a random sequence with varianceσ2=0.49. Choose δ=1.0×10-4. In the proposed algorithm, based on Remark 3.4, select η(0)=0.015,k1=6, θ0=[0,0,0,0,0.1,0.1,0.1,0.1,0]T, and β1(k0)=0, respectively. For comparison, the traditional gradient method is also used to estimate the parameters of the system. In this method, the nonsmooth points of the system are omitted for the gradients of the system do not exist at nonsmooth points. The initialized values of the parameters are the same as those used in the proposed method. The optimizing step is chosen as 0.009.Figure2 illustrates the comparison of the estimated parameter convergence procedures between the proposed method and the traditional gradient method. In Figure 2, blue and solid lines denote the convergence procedures of the parameters estimated by the proposed method, while red and dotted lines show the convergence procedures of the parameters determined by the traditional gradient method. From Figure 2, we note that the parameters of the backlash-like hysteresis converge slower than those of the linear submodels especially the input linear submodel. Moreover, the proposed method has achieved faster convergence than that of the traditional gradient method. It is noticed that the oscillation and sharp jumps happened in the estimation procedure of the traditional gradient approach.Figure 2
The convergence of the estimated parameters (noise-free case): proposed method: blue and solid line; traditional gradient method: red and dotted lines.In the case that the system is affected by random noise, the proposed strategy can still obtain better convergence of parameter estimation. In the simulation with noise, the signal to noise rate (SNR) is equal to 46.5. All the initial values of the parameters are the same as those in the noise-free case.Figure3 shows the comparison of the convergence procedures of the estimated parameters in the case with noise between the proposed method and the traditional gradient approach. Similar to the noise-free case, the blue and solid lines denote the convergence procedures of the parameters estimated by the proposed method, while the red and dotted lines show the convergence procedures of the parameters estimated by the traditional gradient method. Obviously, the proposed method has obtained faster convergent results than the traditional gradient method.Figure 3
The convergence of the estimated parameters (noisy case): proposed method: blue and solid lines; traditional gradient method: red and dotted line.


## 6. Application to anX-Y Moving Positioning Stage

The proposed identification approach is also applied to the modeling of anX-Y moving positioning stage with the architecture shown in Figure 4. In this equipment, the movement of the work platform of each axis is driven by a DC servomotor through a ball-screw-nut mechanism which transforms the rotational shaft movement into linear displacement. The servomotor is controlled by a digital signal processor (TMS320LF-2407A). The displacement of each axis is measured by a linear encoder (RGF2000H125B). The signals of both phase A and phase B encoders are decoded by a quadrature decoding circuit which is based on the decoding chip (Agilent HCTL-2020).Figure 4
TheX-Y moving positioning stage.In this system, the servomotor can be considered as a second-order linear dynamic subsystem. The movement of the work platform is also described by a linear second-order dynamic model. Due to the inherent characteristic, both dead zone and backlash-like hysteresis exist in this system. In order to simplify the identification procedure, the dead zone is compensated by a dead zone inverse model-based compensator. Thus, in the identification, only the effect of backlash-like hysteresis existing in the ball-screw-nut mechanism will be considered. Therefore, the identified system is actually a typical sandwich system with backlash-like hysteresis. In this section, only the identification procedure of axis A will be presented due to the limited space. The corresponding models used to describe the behavior of axis A are shown as follows:(1)
the input linear model (L1):
(6.1)x(k)=-a1x(k-1)-a2x(k-2)+b0u(k-1),(2)
the model of the backlash-like hysteresis:(6.2)v(k)={m1(x(k)-D1),x(k)>v(k-1)m1+D1,x(k)>x(k-1),v(k-1),v(k-1)m2-D2≤x(k)≤v(k-1)m1+D1,m2(x(k)+D2),x(k)<v(k-1)m2-D2,x(k)<x(k-1),(3)
the output linear model (L2):
(6.3)y(k)=-c1y(k-1)-c2y(k-2)+d0v(k)+d1v(k-1),where y(k) is the moving speed of the work platform.Based on the operating requirement, a sequence of square wave plus sinusoidal wave is used to excite the system within the operating range. The corresponding amplitude of the input varies in the range between −1.09 V and 1.05 V, and the sample period is 0.5 ms.In this model, bothb0 and d0 are set to one. The initial values of the other parameters are chosen as η(0)=0.00116, μ=1, θ0=[0,0,0,0,1,1,0.001,0.001,0]T, and β(k0)=0. After 6700 steps, the convergence of the estimation is achieved. Figure 5 illustrates the corresponding procedure of the parameter estimation. It shows that the estimate procedure converges quickly. Figure 6 shows the corresponding mean square error (MSE) of the parameter estimation. We can see that the MSE is decreased sharply in the beginning, at the 180th step, and a local minimum can be found. After that, the algorithm jumps out of the local minimum, and the corresponding MSE gradually converges to a constant of about 0.4.Figure 5
Convergence of parameters of the model of theX-Y moving positioning stage.Figure 6
The MSE plot.Then, the corresponding model validation result is shown in Figure7(a), while Figure 7(b) shows the comparison of the input-output plots between the proposed model and the real data. The maximum relative modeling error is less than 11%. Moreover, it is obvious that the obtained model can accurately approximate the behavior of the X-Y moving positioning stage. Hence, it can be concluded that the proposed identification method is rather promising in engineering application.Model validation result.
(a)
(b)


## 7. Conclusion

In this paper, a recursive gradient-based identification algorithm for the sandwich system with backlash-like hysteresis is proposed. The subgradient is applied to the search of gradient direction at the nonsmooth points of the system. In order to find the proper search direction at the nonsmooth points, the technique of so-called bundle method is utilized. Simulation results have shown that the proposed algorithm has provided us with an option for identification of nonsmooth dynamic systems, and it provides a novel method to identify the more complicated nonsmooth systems. The experimental results ofX-Y stage also show that the proposed method has potential in engineering applications.



---

*Source: 457603-2012-07-15.xml*