- Email: [email protected]
Identification of Hammerstein systems including a nonparametric switched-hysteresis element
Proceedings of the 15th IFAC Symposium on System Identification Saint-Malo, France, July 6-8, 2009
Identification of Hammerstein systems including a nonparametric switched-hysteresis element Y. Rochdi*, F. Giri, F.Z. Chaoui * Corresponding author: GREYC Lab, Université de Caen, France (e-mail: [email protected]) Abstract: Most previous works on Hammerstein system identification focused on the case where the nonlinear element is memoryless. Those which dealt with memory nonlinearities mainly focused on backlash or relay elements. In the present paper, we are considering the case of switched hysteretic nonlinearities. An identification method is developed that yields consistent estimates of the linear subsystem parameters as well as m points on both nonlinearity borders. The latter are allowed to be discontinuous, noninvertible and of unknown structure. 1. INTRODUCTION The problem of identifying Hammerstein systems (Fig 1a) has been given a great deal of interest over the last two decades, see e.g. (Greblicki and Pawlak, 2008) and reference list therein. However, most solutions were developed supposing the system nonlinear element to be memoryless. Furthermore, those solutions that also apply to memory nonlinearities focused mainly on standard backlash- or relayhysteretic elements, (Bai, 2002), (Cerone et al., 2007), (Giri et al., 2008a). The case of switched hysteretic nonlinearities has been recently addressed in (Giri et al., 2008b). The identification scheme proposed there was actually shown to yield consistent estimates of the all unknown system parameters. However, the proposed solution only applies to standard switched hysteresis elements, namely those bordered by two straight lines. In the present paper, we seek the identification of Hammerstein systems that involve arbitrary switched hysteretic nonlinearities. Specifically, the hysteresis ascendant and descendent borders are allowed to be discontinuous, noninvertible and of unknown structure (Fig 2). Such type of hysteresis may be encountered for example, in mechanical structures like suspension system involving leaf spring in heavy vehicles. In this case the hysteresis is created by the dry friction (Coulomb friction) between the leaves. Figure 1.b shows an ideal (kx =cste) model of such leaf spring system and a more realistic, and so complex, one (kx =f(x)). A new identification method is developed that consistently estimates the linear subsystem parameters and m points on both borders of the hysteretic nonlinearity. It involves multiple experiments designed such that, in the corresponding operation conditions: (i) the hysteretic nonlinear element can be assimilated (with no approximation) to a memoryless nonlinearity (memory masking), (ii) parameterizations involving linearly the system unknown parameters can be built-up and based upon to estimate these parameters using least squares estimators, (iii) the resulting signals provide the persistent excitation necessary for the consistency of the above estimators.
978-3-902661-47-0/09/$20.00 © 2009 IFAC
The paper is organized as follows: the identification problem is formally stated in Section 2; the identification scheme design is performed in Sections 3 and 4, a conclusion and reference list end the paper. ξ(t)
Nonlinearity
v(t)
u(t)
F
y(t)
1/A(q-1)
B(q-1)
Linear subsystem
Fig. 1a. Hysteretic Hammerstein System kx
FC fC
fk
f
f = fC + fk f = FC sign( x& ) + k x x
x
x
f
f
f = fC + fk
x
f = FC sign( x& ) + k x ( x ) x
Fig. 1b. Hysteretic leaf spring suspension system 2. IDENTIFICATION PROBLEM STATEMENT 2.1 Class of identified systems We are interested in systems that can be described by the Hammerstein model (Fig. 1a):
342
A(q -1 )y(t) = B(q -1 )u(t) + ξ(t) with
u = Fv
(2.1)
A(q -1 ) = 1 + a 1 q -1 + ... + a n q -n
B(q -1 ) = b1 q -1 + ... + bn q -n
(2.2)
10.3182/20090706-3-FR-2004.0370
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009
{ξ(t)} is a zero-mean stationary sequence of independent random variables. The linear subsystem is supposed to be stable, controllable and of known order n. Controllability is required to get persistent excitation applying a technical lemma in (Giri et al., 2002). In addition to these standard assumptions it is further assumed that B(1)≠0 .
whatever k ≠ 0 . In the next section, we will get a benefit from such model multiplicity imposing a useful condition on the model. u
u
uM
uM
Cd Cd
um Ca
um
v hM
hm
hm
Fig. 2. Nonparametric switched hysteresis The nonlinearity is a switched hysteresis operator characterized by its descendant and ascendant lateral borders C d ( v ) and C a ( v ) . These assume no particular structure unlike in (Bai, 2002), (Cerone et al., 2006), (Giri et al., 2008a-b). The working point (v(t ), u (t )) moves on C d ( v ) and C a ( v ) , but never outside. It jumps from one border to the other each time the increment v(t ) − v(t − 1) changes its sign. Analytically, the nonlinear element F is described as follows: Ca (v(t )) if v(t ) > v(t − 1) u (t ) = u (t − 1) if v(t ) = v(t − 1) C (v(t )) if v(t ) < v(t − 1) d
v
2.2 Identification purpose The purpose is to develop an identification method that accurately determines both the linear subsystem model B(q -1 )/A(q -1 ) and the nonlinear hysteretic operator F . As the latter is nonparametric, the last objective amounts to determine a sufficiently large number of points on each border of F . The two sets of points to be determined are ( V jd , C d (V jd ) ); j = 1,..., m and respectively denoted
{
(2.3)
, C a(V ja
{ } )); j = 1,..., m}. The position of these points will
be précised later.
Assumption 1. There exist known real numbers hm and hM def
def
such that u m < u M with um = C d ( hm ) , u M = C a ( hM ) . Remarks 1. a) Assumption 1 ensures that, when v( t ) spans the working interval [hm , hM ] , the hysteresis output u( t ) takes at least two different values. This is needed for system identifiability. b) The borders Ca (.) and C d (.) assume no particular structure. Therefore, they may be nonparametric, noninvertible and discontinuous. , making it possible to account for rely-type hysteresis. Fig 3 shows a quite general shape of such borders. That is, the present work constitutes a quite significant progress, compared to (Giri et al, 2008b) where the borders are straight lines. c) Finally, it is a well known fact that the system model is not uniquely defined by the representation (2.1). Indeed, if the triplet A( q −1 ), B( q −1 ), F is representative of the system,
)
hM
Fig 3. General shape switched-hysteresis
( V ja
The ascendant border Ca (.) and descendent border C d (.) are only subject to the following assumption:
(
Ca
(
then so is any triplet of the form A( q −1 ),k B( q −1 ), F / k
)
The inaccessibility of the internal sequence u(t) to measurement is a common difficulty to all Hammerstein systems i.e. system identification then must be solely based on the available data, namely the measurements of v( t ) and y( t ) ). A specific difficulty to the case of hysteretic operators
is that the working point (v( t ),u( t )) is likely to move along both hysteresis borders. As a matter of fact, the shape of the trajectory {(v( t ),u( t )) , 0 ≤ t ≤ T0 } (for some finite T0 > 0 ) differs from an experiment to the other. In the present work, periodic signals are privileged, not only for their simplicity, but also because they generate steady-state limit cycles. We seek limit cycles that coincide with one (and only one) border of the hysteresis nonlinearity. Doing so, the latter can be assimilated, during the considered experiment conditions, to a memoryless nonlinearity. 3. IDENTIFICATION OF LINEAR SUBSYSTEM AND SWITCHED-HYSTERESIS DESCENDANT BORDER An identification scheme is determined that determines the linear subsystem parameters ( ai , bi ) as well as m points on the hysteresis descendant border. A key feature of the scheme is the following experiment that makes the working point (v( t ),u( t )) move only on the descendent border.
343
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009
3.1 Descendant Periodic Experiment (DPE) This experiment involves a periodic input signal v d (t ) , with period T = n(m + 3 ) , defined as follows: V d for t = t k + τ j ; j = 1, ..., m v d (t) = j hM otherwise tk = k T
where
( k = 0 ,1, 2 K );
τ 1 = 2n ;
monotonic nature of the sequence (3.1) will prove to be crucial in rendering the signals v d (t ) and u (t ) persistently exciting, making possible a consistent estimation of the above unknown points and parameters.
3.2 System Reparameterization in DPE operation conditions
(3.1)
DPE-model. It was already pointed out that, in DPE operation conditions, the trajectory of v d ( t ) ,u( t ) stays all time on Dlc , occupying there the specific positions
(
τ j = ( 2 + j)n
( j = 2 , ..., m ). The V ’s are any real numbers such that d j
(hM
hm = V1d < V2d < ... < Vmd < hM . It is readily seen that the input
{
}
signal v d ( t ) is pulse type: in each period [t k ,t k +1 ) , it equals hM all time except at the instants tk + τ j ( j = 1, 2 , ..., m )
where it takes different values, namely V1d , V2d ,..., Vmd , which all are smaller than hM . It follows from (2.3) that the
(
)
working point v d (t ),u (t )
stays all time on the limit cycle,
denoted Dlc , including the single point ( hM ,u M ) and the descendant branch
{(v ,C ( v )), h d
m
}
≤ v < hM . It is clear that
the domain Dlc can also be seen as a simple curve that uniquely defines a standard function, denoted F d , such that Dlc = v , F d ( v ) , hm ≤ v ≤ hM . Specifically:
{(
)
}
C ( v ) if hm ≤ v < hM F (v) = d u M if v = hM d
(3.2)
During the present experiment, the working point v d ( t ) ,u( t ) never gets to the ‘ascending’ branch {(v ,Ca ( v )) , hm ≤ v < hM }; it stays all time on Dlc . Furthermore, it never goes directly from one position V jd ,Cd ( V jd ) to the next one V jd+1 ,Cd ( V jd+1 ) ; it necessarily
(
)
(
)
(
)
d
(V
d d j +1 , C d (V j +1 )
(V
d j
, Cd (V jd )
)
, F d (V jd )
)
( j = 1,..., m) . Then, the
(
)
A(q -1 )y(t) = B(q -1 )F d (v(t)) + ξ(t) with v(t) = v d (t) (3.3)
where F d is defined by (3.2). Note that the substitution of F d (.) to F (.) involves no error as long as the applied input sequence is v d (t) (defined by (3.1)).
(
)
The fact that the couple v d (t ) , F d (v d (t )) only occupies the
(hM
positions
,u M ) and
(
V jd
d
, F (V jd )
)
( j = 1,..., m) ,
makes possible an exact polynomial representation of the (nonparametric) nonlinearity F d . Before performing such polynomial representation, we will first proceed with the centering of F d so that the obtained centered nonlinearity, ~ ~ denoted Fd , satisfies Fd ( 0 ) = 0 .
{
}
DPE-Model Centering. Let y0d ( t ) denote the response of (2.1) to the following step input: def h v(t) = v0d ( t ) = m hM
if t = 0 if t ≥ 1
(3.4)
Indeed, one gets from (2.3) (or Fig. 2) that, for all t ≥ 1 :
Cd
uM
(
d
d j
quantity F v (t) may be substituted to F v d (t) in the initial model (2.1). Doing so, one gets a new model, referred to DPE-model, defined by:
transits via the ‘satellite’ position ( hM ,u M ) , see Fig 4. u
(V
,u M ) and
)
def
Satellite
)
(
)
u0d ( t ) = F v0d ( t ) = u M whatever u0d ( 0 ) . As the point (hM ,u M ) belongs to the domain Dlc , one also has def
(
)
u0d ( t ) = F d v0d ( t ) = u M . Therefore, the couple of signals
)
{v
Ca
Vmd
hM v
{
}
A(q -1 )y 0d (t ) = B(q -1 )u 0d (t ) + ξ 0d (t )
Fig. 4. Movement of the working point (v(t), u(t)) on the limit cycle Dlc during the descendent periodic experiment The data collected in the current experiment will be used to determine the points V jd , F d V jd ( j = 1,2...,m ), as well as
(
}
d 0
( t ) and y ( t ) can be related by the general model (2.1) as well as by the DPE-model (3.3). Thus, one has for all t ≥ 1:
um
V1d = hm V2d
d 0
def
(
)
(3.5)
where u 0d (t ) = F d v0d (t ) = u M and ξ 0d ( t ) denotes the realization of ξ ( t ) during the present experiment. Time-
( ))
the linear subsystem parameters ( ai , bi ). The pulse and non344
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009
averaging the first equality in (3.3), over the interval 1≤ t≤ L, yields1: A( 1 ) y 0d ( L ) = B( 1 ) u M + ξ 0d ( L )
(3.6)
As ξ 0d ( t ) is zero mean and ergodic, one has ξ 0d ( L ) → 0 as L → ∞ . Then, it follows from (3.6) that: A( 1 ) y 0d = B( 1 ) u M
where y
(3.7)
is the limit of y ( L ) when L → ∞ . Practically,
d 0
d 0
y is obtained from a sufficient data sample {y(t), 1 ≤ t ≤ L} . d 0
~ If the original F d does not satisfy (3.11c), consider the ~ model A( q −1 ),k B( q −1 ), F d / k with k = −(um − u M ) . Indeed, the latter does satisfy (3.11c), due to (3.11a-b). Now, recall that the initial identification purpose is to determine, on one hand, the coefficients of the polynomial pair A( q −1 ) , B( q −1 ) and, on the other hand, the m points
(
{ {(V
(
)
(
)
d j
)
}
of
Fd
the
graphical
{(W
d j
)
}
def ~ , F d (W jd ) ; j = 1,2..., m with W jd = V jd − hM (3.12)
To
this
end,
we
v( t ) = v d ( t )
let
so
that
def
(3.8)
v~( t ) = v~ d ( t ) = v d ( t ) − hM is in turn a pulse-type periodic
v~( t ) = v( t ) − hM , ~y ( t ) = y( t ) − y 0d ,
(3.9a)
signal with the same period as v d ( t ) . Specifically, one has from (3.2a) and (3.12) that, for t k ≤ t < t k +1 (with tk = k T and k = 0 ,1,2K ):
def ~ F d ( w ) = F d ( w + hM ) − u M
(3.9b)
For convenience, let us introduce the notations: def
def
W d for t = t k + τ j ; j = 1,2...,m v~ d ( t ) = j 0 otherwise
def ~ It is readily seen that F d (v~(t)) = F d (v(t)) − u M . Then, (3.8)
can be given the compact form, called centered DPE- model: ~ A( q ) ~ y ( t ) = B( q −1 ) F d (v~( t )) + ξ ( t ) −1
~ ~ F d (0) = 0 , F d (W1d ) = u m − u M ≠ 0 def
W1d = V1d − hM = hm − hM
(3.11a)
where the fact that um − u M ≠ 0 is ensured by Assumption 1. Properties (3.11a) will prove to be crucial in achieving an important persistent excitation property (see Subsection 3.3). For now, let us use (3.11a) to get off model multiplicity pointed out in Remarks 1 (Part c). Accordingly, any triplet of ~ the form A( q −1 ),k B( q −1 ), F d / k is also representative of the model (3.10), whatever k ≠ 0 . In the latter, the focus will be made on the unique model that satisfies:
(
}
)
(
}
)
(
)
~ Polynomial representation of F d v~ d ( t ) . Let P d ( w ) denotes the unique mth degree polynomial interpolating the ~ points (0,0) and W jd , F d ( W jd ) ; ( j = 1,2 ,.., m ) . Such
)
polynomial can simply be described using the Vandermonde formula: m
P d ( w) = w ∑ c dj Pjd ( w)
(3.14a)
j =1
)
~ F d (W1d ) = − 1
{
values in the set 0 ,W jd ; j = 1, 2 , ..., m and, consequently, ~ the trajectory v~ d ( t ) , F d ( v~ d ( t )) ; t ≥ 0 consists of m+1 points: the m points of the set in (3.12) and the origin (0,0). ~ Therefore, the function F d v d ( t ) can be given an exact polynomial representation using a polynomial of degree m.
(
(3.11b)
(3.13)
where the τ j ’s are as in (3.1). Clearly, v~ d ( t ) takes its
{(
(3.10)
Equation (3.10) defines a Hammerstein model involving the ~ nonlinearity F d (.) which is analytically defined in (3.9b). Using (3.2) (or Fig. 4), one readily gets from (3.9b), that:
1
}
, F ( V ) , j = 1,2 ,..., m d
characteristic. It follows from (3.9b) that the last requirement ~ amounts to determining m points of the F d characteristic, namely:
Subtracting both sides of (3.7) from corresponding sides of DPE-model (3.3), gives: A(q -1 ) y(t) − y 0d = B(q -1 ) F d (v(t)) − u M + ξ(t) (for all t ≥ 0 )
d j
)
c dj
=
~ F d (W jd ) W jd
def m
; Pjd ( w) = Π
w − Wi d
i =1 W d j i≠ j
(3.14b)
− Wi d
It can be easily seen that:
(3.11c)
1 Pjd ( Wi d ) = δ ij = 0
Throughout the paper, x ( L ) denotes the arithmetic average of a L −1 sequence x( t ) over the time interval [0 L] i.e. x ( L ) = 1 ∑ x( t ) . L t =0
In the case of stationary ergodic random sequences, one has E ( x( t )) = lim x ( L ) , for all t . L→ ∞
345
if
i= j
if
i≠ j
(3.15)
~ P d ( 0 ) = 0; P d (Wjd ) = F d (Wjd ) ; j = 1,...,m
Then,
{0 ,W
d j
as
v~ d ( t )
}
; j = 1,2 ,...,m ,
only
takes it
values
(3.16) in
follows
the
set that
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009
(
)
(
)
(
)
~ F d v~ d ( t ) = P d v~ d ( t ) , for all t. Therefore, P d v~ d ( t ) can ~ be substituted to F d v~ d ( t ) in the D-Centered model (3.8). Doing so, one gets, (for all t ≥ 0 ):
(
)
A( q −1 ) ~ y ( t ) = B( q −1 )P d ( v~ d ( t )) + ξ ( t )
DPE Regression Form. Using (3.14a), one gets from (3.17): A(q ) ~ y (t ) = −1
∑∑
i =1 j =1
µ ijd
(
together with (3.9b), (3.11b) and (3.7) implies: F d (V jd ) = c dj W jd +
(3.17)
The main advantage of (3.17), over (3.10), is that it involves a parametric nonlinearity making possible the construction of a regression form. It is worthy pointing out that the passage from (3.10) to (3.17) generates no error, as long as the used input is the signal v~ d ( t ) which satisfies (3.13).
n m
~ Now, (3.14) immediately gives F d (W jd ) = c dj W jd which,
Estimation of the ai ’s and µijd ’s. Based on the regression
(3.19a-c), the unknown parameter vector Θ d can be estimated simply using the standard least squares estimator: 1 N
Θˆ d ( N ) =
)
d d µ11 ... µ nd1 ... µ1dm ... µ nm
]
] (3.19b)
~ c1d = F d ( W1d ) / W1d = −1 / W1d
T
(3.19c) Proposition 3.1 (Consistency of estimator (3.21). Consider the model (3.17) and its undisturbed version (3.22), both excited by the signal v~ d ( t ) defined by (3.13). Then one has the following properties:
1) The state vector Z d ( t ) , generated by (3.22a-b), possesses the strong persistent excitation property i.e. there exists a real λ >0 such that, for all k: 4 n −1
(3.20a)
∑Z
(i=1, ..., n)
(3.20b)
lim
N →∞
bi and summing over i = 1,K ,n :
∑ i =1
bi2
(j=2, ..., m)
(3.20c)
( t k + i )Z d ( t k + i )T ≥ λ I
(I is the identity matrix)
i =0
The c dj ’s are in turn obtained simply multiplying (3.18b) by
d ij
d
2) The state vector Φ d ( t ) generated by (3.17) (or equally (3.10)) possesses the weak persistent excitation property i.e. there exists a real β>0, such that:
W1d ’s). Then, using (3.18b) one gets the bi ’s:
i
(3.22a)
]
That is, the coefficient c1d is perfectly known (because so is
i =1 n
)
v~ d ( t − 1 )Pmd ( v~ d ( t − 1 )) ...v~ d ( t − n )Pmd ( v~ d ( t − n )) (3.22b)
This question is answered making use of available properties. First, it readily follows from (3.11c) and (3.14b):
c =
(3.21)
v~ d ( t − 1 )P1d ( v~ d ( t − 1 )) ... v~ d ( t − n )P1d ( v~ d ( t − n )) ...
we obtain estimates of ( bi , c dj ) from those of µijd ?
∑b µ
( i )~ y ( i )T
[
As Θ d comes in linearly, equation (3.19c) turns out to be an adequate parameterization to get estimates of the parameters ai and µijd , using a least squares type algorithm. But, how can
d j
i =1
d
Z d ( t )T = − z d ( t − 1 )... − z d ( t − n )
v~ d (t − 1) Pmd (v~ d (t − 1))...v~ d (t − n) Pmd (v~ d (t − n))
n
N
∑Φ
Introduce the undisturbed regression vector:
(3.19a)
v~ d (t − 1) P1d (v~ d (t − 1)) ... v~ d (t − n) P1d (v~ d (t − n)) ...
bi = µ id1 / c1d = − µ id1 W1d
1 ( i )Φ d ( i )T N
(
Φ d (t ) T = [ − ~y (t − 1)... − ~y (t − n)
[
i =1
−1
d
A( q −1 ) z d ( t ) = B( q −1 )P d v~ d ( t )
Equation (3.18a) is given the linear regression form:
Θ d = a1 ... a n
N
∑Φ
This estimator will now be shown to be consistent. Let us consider the undisturbed version of model (3.17):
(3.18b)
~ y (t) = Φ d (t)T Θ d + ξ ( t )
(3.20d)
3.3 Parameter estimation
v~ d (t − i ) Pjd v~ d (t − i ) + ξ (t ) (3.18a)
µ ijd = bi c dj (i = 1,..., n; j = 1,..., m )
A(1) d y0 B (1)
1 N
N
∑Φ
d
( t )Φ d ( t )T > β I
(w.p.1)
i =1
3) The estimator (3.21) is consistent i.e. Θˆ d ( N ) → Θ d as N → ∞ (w.p. 1) Proof. The proof is omitted due to the limitation on the paper length. Its main ingredients can be found in (Chaoui et al., 2006).
346
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009 tk = k T
( k = 0 ,1,2K );
τ1 = 2n ;
τ j = ( 2 + j)n
Estimation of the bi ’s and the c j ’s. Relations (3.20b-c)
where
suggest the following estimators:
( j = 2 , ..., m ). The resulting limit cycle, denoted Alc , is shown in (Fig. 5). The rest of identification scheme development follows closely Section 3. The detailed design and analysis is omitted due to space limitation.
bˆi ( N ) = − µˆ i 1 ( N )W1d
(i=1, ..., n)
(3.24a)
(j=2, ..., m)
(3.24b)
n
∑ bˆ ( N )µˆ ( N ) i
cˆ dj ( N ) =
u
ij
i =1
n
∑ bˆ ( N ) 2 i
Cd
uM
i =1
(V
Proposition 3.2. The estimators (3.24a-b) are consistent. Proof. Comparing (3.24a) and (3.20b) yields ˆb ( N ) − b = (µ − µˆ ( N ))W d . The consistency of (3.24a) is i i i1 i1 1 then a direct consequence of the consistency of Θˆ d ( N )
(because this in particular ensures that µi 1 − µˆ i 1 ( N ) converges in probability to zero as N → ∞ ). The consistency of (3.24b) is proved similarly using the fact that both µˆ i 1 ( N ) and bˆi ( N ) are consistent estimators Estimation
of
the
points
(V
d j
)
Aˆ (q −1 , N ), Bˆ (q −1 , N )
A(q −1 ), B(q −1 )
are
the
estimates
by the previously obtained parameter estimates aˆ i ( N ), bˆi ( N ) . Comparing (3.25) and (3.20d) implies: Fˆ jd ( V jd , N ) − F d ( V jd ) =
(
)
hm
ˆ ( 1, N ) d A A( 1 ) d . y0 ( N ) − y0 B( 1 ) Bˆ( 1, N )
The consistency of (3.25) then follows from the consistency of cˆ dj ( N ), aˆ i ( N ), bˆi ( N ) and y 0d ( N ) . 4. IDENTIFICATION OF THE LINEAR SUBSYSTEM AND SWITCHED-HYSTERESIS ASCENDANT BORDER An identification scheme is developed that determines the linear subsystem parameters ( ai , bi ) as well as m points on the hysteresis ascendant border. There is a strong symmetry between the structure of this identification scheme and the one described in section 3. A key ingredient of the new scheme is an ascendant periodic experiment (APE) involving the periodic input signal v a ( t ) , with period T = n(m + 3 ) , defined as follows: V a for t = tk + τ j ; j = 1, ..., m v a(t) = j hm otherwise
)
a V2a V1 = hvM
V ma
Fig. 5. Movement of the working point on the limit cycle Alc when using the periodic input signal (4.1) 5. CONCLUSION
REFERENCES
of
induced
W jd cˆ dj ( N ) − c dj +
,C a (V jd+1 )
)
Hammerstein systems identification is considered is presence of memory nonlinearities of the switched-hysteresis type. The nonlinearity borders may be discontinuous, noninvertible and of unknown structure.
Aˆ (1, N ) d Fˆ d (V jd , N ) = cˆ dj ( N )W jd + y 0 ( N ) (3.25) Bˆ (1, N )
where
d j +1
,C a(V ja )
Ca
um
, F d (V jd ) , j = 2,..., m .
Equation (3.20d) suggests for F d ( V jd ) the estimator:
(V
Satellite
a j
Bai, E.W. (2002).‘Identification of linear systems with hard input nonlinearities of known structure’. Automatica, vol 38, pp. 853-860. Cerone, V., Regruto, D. (2007) ‘Bounding the parameters of linear systems with input backlash’. IEEE Trans. Automatic Control, vol 52(3), pp.531–536. Chaoui, F.Z., Giri, F., Rochdi, Y., Haloua, M., Naitali, A.(2005). ‘System identification based on Hammerstein model’. International Journal of Control, vol. 78, N° 6, pp. 430-442. Giri, F., Chaoui, F.Z., Rochdi, Y. (2002). ‘Interval excitation through impulse sequences. A technical lemma’. Automatica, vol 38, pp. 457-465. Giri, F., Rochdi, Y., Chaoui, F.Z., Brouri. A. (2008a). ‘Identification of Hammerstein systems in presence of hysteresis-backlash and hysteresis-relay nonlinearities’. Automatica, vol. 44, pp 767-775. Giri, F., Rochdi, Y., Elayan, E., Brouri, A., Chaoui, F.Z. (2008b) ‘Hammerstein systems identification in presence of hysteresis-backlash nonlinearity’. IFAC World Congress, Seoul, pp. 7859-7864. Greblicki, W., M. Pawlak (2008). ‘Nonparametric System Identification’. Cambridge University Press. 2008. ISBN-10: 0521868041.
(4.1)
347
Identification of Hammerstein systems including a nonparametric switched-hysteresis element Y. Rochdi*, F. Giri, F.Z. Chaoui * Corresponding author: GREYC Lab, Université de Caen, France (e-mail: [email protected]) Abstract: Most previous works on Hammerstein system identification focused on the case where the nonlinear element is memoryless. Those which dealt with memory nonlinearities mainly focused on backlash or relay elements. In the present paper, we are considering the case of switched hysteretic nonlinearities. An identification method is developed that yields consistent estimates of the linear subsystem parameters as well as m points on both nonlinearity borders. The latter are allowed to be discontinuous, noninvertible and of unknown structure. 1. INTRODUCTION The problem of identifying Hammerstein systems (Fig 1a) has been given a great deal of interest over the last two decades, see e.g. (Greblicki and Pawlak, 2008) and reference list therein. However, most solutions were developed supposing the system nonlinear element to be memoryless. Furthermore, those solutions that also apply to memory nonlinearities focused mainly on standard backlash- or relayhysteretic elements, (Bai, 2002), (Cerone et al., 2007), (Giri et al., 2008a). The case of switched hysteretic nonlinearities has been recently addressed in (Giri et al., 2008b). The identification scheme proposed there was actually shown to yield consistent estimates of the all unknown system parameters. However, the proposed solution only applies to standard switched hysteresis elements, namely those bordered by two straight lines. In the present paper, we seek the identification of Hammerstein systems that involve arbitrary switched hysteretic nonlinearities. Specifically, the hysteresis ascendant and descendent borders are allowed to be discontinuous, noninvertible and of unknown structure (Fig 2). Such type of hysteresis may be encountered for example, in mechanical structures like suspension system involving leaf spring in heavy vehicles. In this case the hysteresis is created by the dry friction (Coulomb friction) between the leaves. Figure 1.b shows an ideal (kx =cste) model of such leaf spring system and a more realistic, and so complex, one (kx =f(x)). A new identification method is developed that consistently estimates the linear subsystem parameters and m points on both borders of the hysteretic nonlinearity. It involves multiple experiments designed such that, in the corresponding operation conditions: (i) the hysteretic nonlinear element can be assimilated (with no approximation) to a memoryless nonlinearity (memory masking), (ii) parameterizations involving linearly the system unknown parameters can be built-up and based upon to estimate these parameters using least squares estimators, (iii) the resulting signals provide the persistent excitation necessary for the consistency of the above estimators.
978-3-902661-47-0/09/$20.00 © 2009 IFAC
The paper is organized as follows: the identification problem is formally stated in Section 2; the identification scheme design is performed in Sections 3 and 4, a conclusion and reference list end the paper. ξ(t)
Nonlinearity
v(t)
u(t)
F
y(t)
1/A(q-1)
B(q-1)
Linear subsystem
Fig. 1a. Hysteretic Hammerstein System kx
FC fC
fk
f
f = fC + fk f = FC sign( x& ) + k x x
x
x
f
f
f = fC + fk
x
f = FC sign( x& ) + k x ( x ) x
Fig. 1b. Hysteretic leaf spring suspension system 2. IDENTIFICATION PROBLEM STATEMENT 2.1 Class of identified systems We are interested in systems that can be described by the Hammerstein model (Fig. 1a):
342
A(q -1 )y(t) = B(q -1 )u(t) + ξ(t) with
u = Fv
(2.1)
A(q -1 ) = 1 + a 1 q -1 + ... + a n q -n
B(q -1 ) = b1 q -1 + ... + bn q -n
(2.2)
10.3182/20090706-3-FR-2004.0370
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009
{ξ(t)} is a zero-mean stationary sequence of independent random variables. The linear subsystem is supposed to be stable, controllable and of known order n. Controllability is required to get persistent excitation applying a technical lemma in (Giri et al., 2002). In addition to these standard assumptions it is further assumed that B(1)≠0 .
whatever k ≠ 0 . In the next section, we will get a benefit from such model multiplicity imposing a useful condition on the model. u
u
uM
uM
Cd Cd
um Ca
um
v hM
hm
hm
Fig. 2. Nonparametric switched hysteresis The nonlinearity is a switched hysteresis operator characterized by its descendant and ascendant lateral borders C d ( v ) and C a ( v ) . These assume no particular structure unlike in (Bai, 2002), (Cerone et al., 2006), (Giri et al., 2008a-b). The working point (v(t ), u (t )) moves on C d ( v ) and C a ( v ) , but never outside. It jumps from one border to the other each time the increment v(t ) − v(t − 1) changes its sign. Analytically, the nonlinear element F is described as follows: Ca (v(t )) if v(t ) > v(t − 1) u (t ) = u (t − 1) if v(t ) = v(t − 1) C (v(t )) if v(t ) < v(t − 1) d
v
2.2 Identification purpose The purpose is to develop an identification method that accurately determines both the linear subsystem model B(q -1 )/A(q -1 ) and the nonlinear hysteretic operator F . As the latter is nonparametric, the last objective amounts to determine a sufficiently large number of points on each border of F . The two sets of points to be determined are ( V jd , C d (V jd ) ); j = 1,..., m and respectively denoted
{
(2.3)
, C a(V ja
{ } )); j = 1,..., m}. The position of these points will
be précised later.
Assumption 1. There exist known real numbers hm and hM def
def
such that u m < u M with um = C d ( hm ) , u M = C a ( hM ) . Remarks 1. a) Assumption 1 ensures that, when v( t ) spans the working interval [hm , hM ] , the hysteresis output u( t ) takes at least two different values. This is needed for system identifiability. b) The borders Ca (.) and C d (.) assume no particular structure. Therefore, they may be nonparametric, noninvertible and discontinuous. , making it possible to account for rely-type hysteresis. Fig 3 shows a quite general shape of such borders. That is, the present work constitutes a quite significant progress, compared to (Giri et al, 2008b) where the borders are straight lines. c) Finally, it is a well known fact that the system model is not uniquely defined by the representation (2.1). Indeed, if the triplet A( q −1 ), B( q −1 ), F is representative of the system,
)
hM
Fig 3. General shape switched-hysteresis
( V ja
The ascendant border Ca (.) and descendent border C d (.) are only subject to the following assumption:
(
Ca
(
then so is any triplet of the form A( q −1 ),k B( q −1 ), F / k
)
The inaccessibility of the internal sequence u(t) to measurement is a common difficulty to all Hammerstein systems i.e. system identification then must be solely based on the available data, namely the measurements of v( t ) and y( t ) ). A specific difficulty to the case of hysteretic operators
is that the working point (v( t ),u( t )) is likely to move along both hysteresis borders. As a matter of fact, the shape of the trajectory {(v( t ),u( t )) , 0 ≤ t ≤ T0 } (for some finite T0 > 0 ) differs from an experiment to the other. In the present work, periodic signals are privileged, not only for their simplicity, but also because they generate steady-state limit cycles. We seek limit cycles that coincide with one (and only one) border of the hysteresis nonlinearity. Doing so, the latter can be assimilated, during the considered experiment conditions, to a memoryless nonlinearity. 3. IDENTIFICATION OF LINEAR SUBSYSTEM AND SWITCHED-HYSTERESIS DESCENDANT BORDER An identification scheme is determined that determines the linear subsystem parameters ( ai , bi ) as well as m points on the hysteresis descendant border. A key feature of the scheme is the following experiment that makes the working point (v( t ),u( t )) move only on the descendent border.
343
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009
3.1 Descendant Periodic Experiment (DPE) This experiment involves a periodic input signal v d (t ) , with period T = n(m + 3 ) , defined as follows: V d for t = t k + τ j ; j = 1, ..., m v d (t) = j hM otherwise tk = k T
where
( k = 0 ,1, 2 K );
τ 1 = 2n ;
monotonic nature of the sequence (3.1) will prove to be crucial in rendering the signals v d (t ) and u (t ) persistently exciting, making possible a consistent estimation of the above unknown points and parameters.
3.2 System Reparameterization in DPE operation conditions
(3.1)
DPE-model. It was already pointed out that, in DPE operation conditions, the trajectory of v d ( t ) ,u( t ) stays all time on Dlc , occupying there the specific positions
(
τ j = ( 2 + j)n
( j = 2 , ..., m ). The V ’s are any real numbers such that d j
(hM
hm = V1d < V2d < ... < Vmd < hM . It is readily seen that the input
{
}
signal v d ( t ) is pulse type: in each period [t k ,t k +1 ) , it equals hM all time except at the instants tk + τ j ( j = 1, 2 , ..., m )
where it takes different values, namely V1d , V2d ,..., Vmd , which all are smaller than hM . It follows from (2.3) that the
(
)
working point v d (t ),u (t )
stays all time on the limit cycle,
denoted Dlc , including the single point ( hM ,u M ) and the descendant branch
{(v ,C ( v )), h d
m
}
≤ v < hM . It is clear that
the domain Dlc can also be seen as a simple curve that uniquely defines a standard function, denoted F d , such that Dlc = v , F d ( v ) , hm ≤ v ≤ hM . Specifically:
{(
)
}
C ( v ) if hm ≤ v < hM F (v) = d u M if v = hM d
(3.2)
During the present experiment, the working point v d ( t ) ,u( t ) never gets to the ‘ascending’ branch {(v ,Ca ( v )) , hm ≤ v < hM }; it stays all time on Dlc . Furthermore, it never goes directly from one position V jd ,Cd ( V jd ) to the next one V jd+1 ,Cd ( V jd+1 ) ; it necessarily
(
)
(
)
(
)
d
(V
d d j +1 , C d (V j +1 )
(V
d j
, Cd (V jd )
)
, F d (V jd )
)
( j = 1,..., m) . Then, the
(
)
A(q -1 )y(t) = B(q -1 )F d (v(t)) + ξ(t) with v(t) = v d (t) (3.3)
where F d is defined by (3.2). Note that the substitution of F d (.) to F (.) involves no error as long as the applied input sequence is v d (t) (defined by (3.1)).
(
)
The fact that the couple v d (t ) , F d (v d (t )) only occupies the
(hM
positions
,u M ) and
(
V jd
d
, F (V jd )
)
( j = 1,..., m) ,
makes possible an exact polynomial representation of the (nonparametric) nonlinearity F d . Before performing such polynomial representation, we will first proceed with the centering of F d so that the obtained centered nonlinearity, ~ ~ denoted Fd , satisfies Fd ( 0 ) = 0 .
{
}
DPE-Model Centering. Let y0d ( t ) denote the response of (2.1) to the following step input: def h v(t) = v0d ( t ) = m hM
if t = 0 if t ≥ 1
(3.4)
Indeed, one gets from (2.3) (or Fig. 2) that, for all t ≥ 1 :
Cd
uM
(
d
d j
quantity F v (t) may be substituted to F v d (t) in the initial model (2.1). Doing so, one gets a new model, referred to DPE-model, defined by:
transits via the ‘satellite’ position ( hM ,u M ) , see Fig 4. u
(V
,u M ) and
)
def
Satellite
)
(
)
u0d ( t ) = F v0d ( t ) = u M whatever u0d ( 0 ) . As the point (hM ,u M ) belongs to the domain Dlc , one also has def
(
)
u0d ( t ) = F d v0d ( t ) = u M . Therefore, the couple of signals
)
{v
Ca
Vmd
hM v
{
}
A(q -1 )y 0d (t ) = B(q -1 )u 0d (t ) + ξ 0d (t )
Fig. 4. Movement of the working point (v(t), u(t)) on the limit cycle Dlc during the descendent periodic experiment The data collected in the current experiment will be used to determine the points V jd , F d V jd ( j = 1,2...,m ), as well as
(
}
d 0
( t ) and y ( t ) can be related by the general model (2.1) as well as by the DPE-model (3.3). Thus, one has for all t ≥ 1:
um
V1d = hm V2d
d 0
def
(
)
(3.5)
where u 0d (t ) = F d v0d (t ) = u M and ξ 0d ( t ) denotes the realization of ξ ( t ) during the present experiment. Time-
( ))
the linear subsystem parameters ( ai , bi ). The pulse and non344
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009
averaging the first equality in (3.3), over the interval 1≤ t≤ L, yields1: A( 1 ) y 0d ( L ) = B( 1 ) u M + ξ 0d ( L )
(3.6)
As ξ 0d ( t ) is zero mean and ergodic, one has ξ 0d ( L ) → 0 as L → ∞ . Then, it follows from (3.6) that: A( 1 ) y 0d = B( 1 ) u M
where y
(3.7)
is the limit of y ( L ) when L → ∞ . Practically,
d 0
d 0
y is obtained from a sufficient data sample {y(t), 1 ≤ t ≤ L} . d 0
~ If the original F d does not satisfy (3.11c), consider the ~ model A( q −1 ),k B( q −1 ), F d / k with k = −(um − u M ) . Indeed, the latter does satisfy (3.11c), due to (3.11a-b). Now, recall that the initial identification purpose is to determine, on one hand, the coefficients of the polynomial pair A( q −1 ) , B( q −1 ) and, on the other hand, the m points
(
{ {(V
(
)
(
)
d j
)
}
of
Fd
the
graphical
{(W
d j
)
}
def ~ , F d (W jd ) ; j = 1,2..., m with W jd = V jd − hM (3.12)
To
this
end,
we
v( t ) = v d ( t )
let
so
that
def
(3.8)
v~( t ) = v~ d ( t ) = v d ( t ) − hM is in turn a pulse-type periodic
v~( t ) = v( t ) − hM , ~y ( t ) = y( t ) − y 0d ,
(3.9a)
signal with the same period as v d ( t ) . Specifically, one has from (3.2a) and (3.12) that, for t k ≤ t < t k +1 (with tk = k T and k = 0 ,1,2K ):
def ~ F d ( w ) = F d ( w + hM ) − u M
(3.9b)
For convenience, let us introduce the notations: def
def
W d for t = t k + τ j ; j = 1,2...,m v~ d ( t ) = j 0 otherwise
def ~ It is readily seen that F d (v~(t)) = F d (v(t)) − u M . Then, (3.8)
can be given the compact form, called centered DPE- model: ~ A( q ) ~ y ( t ) = B( q −1 ) F d (v~( t )) + ξ ( t ) −1
~ ~ F d (0) = 0 , F d (W1d ) = u m − u M ≠ 0 def
W1d = V1d − hM = hm − hM
(3.11a)
where the fact that um − u M ≠ 0 is ensured by Assumption 1. Properties (3.11a) will prove to be crucial in achieving an important persistent excitation property (see Subsection 3.3). For now, let us use (3.11a) to get off model multiplicity pointed out in Remarks 1 (Part c). Accordingly, any triplet of ~ the form A( q −1 ),k B( q −1 ), F d / k is also representative of the model (3.10), whatever k ≠ 0 . In the latter, the focus will be made on the unique model that satisfies:
(
}
)
(
}
)
(
)
~ Polynomial representation of F d v~ d ( t ) . Let P d ( w ) denotes the unique mth degree polynomial interpolating the ~ points (0,0) and W jd , F d ( W jd ) ; ( j = 1,2 ,.., m ) . Such
)
polynomial can simply be described using the Vandermonde formula: m
P d ( w) = w ∑ c dj Pjd ( w)
(3.14a)
j =1
)
~ F d (W1d ) = − 1
{
values in the set 0 ,W jd ; j = 1, 2 , ..., m and, consequently, ~ the trajectory v~ d ( t ) , F d ( v~ d ( t )) ; t ≥ 0 consists of m+1 points: the m points of the set in (3.12) and the origin (0,0). ~ Therefore, the function F d v d ( t ) can be given an exact polynomial representation using a polynomial of degree m.
(
(3.11b)
(3.13)
where the τ j ’s are as in (3.1). Clearly, v~ d ( t ) takes its
{(
(3.10)
Equation (3.10) defines a Hammerstein model involving the ~ nonlinearity F d (.) which is analytically defined in (3.9b). Using (3.2) (or Fig. 4), one readily gets from (3.9b), that:
1
}
, F ( V ) , j = 1,2 ,..., m d
characteristic. It follows from (3.9b) that the last requirement ~ amounts to determining m points of the F d characteristic, namely:
Subtracting both sides of (3.7) from corresponding sides of DPE-model (3.3), gives: A(q -1 ) y(t) − y 0d = B(q -1 ) F d (v(t)) − u M + ξ(t) (for all t ≥ 0 )
d j
)
c dj
=
~ F d (W jd ) W jd
def m
; Pjd ( w) = Π
w − Wi d
i =1 W d j i≠ j
(3.14b)
− Wi d
It can be easily seen that:
(3.11c)
1 Pjd ( Wi d ) = δ ij = 0
Throughout the paper, x ( L ) denotes the arithmetic average of a L −1 sequence x( t ) over the time interval [0 L] i.e. x ( L ) = 1 ∑ x( t ) . L t =0
In the case of stationary ergodic random sequences, one has E ( x( t )) = lim x ( L ) , for all t . L→ ∞
345
if
i= j
if
i≠ j
(3.15)
~ P d ( 0 ) = 0; P d (Wjd ) = F d (Wjd ) ; j = 1,...,m
Then,
{0 ,W
d j
as
v~ d ( t )
}
; j = 1,2 ,...,m ,
only
takes it
values
(3.16) in
follows
the
set that
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009
(
)
(
)
(
)
~ F d v~ d ( t ) = P d v~ d ( t ) , for all t. Therefore, P d v~ d ( t ) can ~ be substituted to F d v~ d ( t ) in the D-Centered model (3.8). Doing so, one gets, (for all t ≥ 0 ):
(
)
A( q −1 ) ~ y ( t ) = B( q −1 )P d ( v~ d ( t )) + ξ ( t )
DPE Regression Form. Using (3.14a), one gets from (3.17): A(q ) ~ y (t ) = −1
∑∑
i =1 j =1
µ ijd
(
together with (3.9b), (3.11b) and (3.7) implies: F d (V jd ) = c dj W jd +
(3.17)
The main advantage of (3.17), over (3.10), is that it involves a parametric nonlinearity making possible the construction of a regression form. It is worthy pointing out that the passage from (3.10) to (3.17) generates no error, as long as the used input is the signal v~ d ( t ) which satisfies (3.13).
n m
~ Now, (3.14) immediately gives F d (W jd ) = c dj W jd which,
Estimation of the ai ’s and µijd ’s. Based on the regression
(3.19a-c), the unknown parameter vector Θ d can be estimated simply using the standard least squares estimator: 1 N
Θˆ d ( N ) =
)
d d µ11 ... µ nd1 ... µ1dm ... µ nm
]
] (3.19b)
~ c1d = F d ( W1d ) / W1d = −1 / W1d
T
(3.19c) Proposition 3.1 (Consistency of estimator (3.21). Consider the model (3.17) and its undisturbed version (3.22), both excited by the signal v~ d ( t ) defined by (3.13). Then one has the following properties:
1) The state vector Z d ( t ) , generated by (3.22a-b), possesses the strong persistent excitation property i.e. there exists a real λ >0 such that, for all k: 4 n −1
(3.20a)
∑Z
(i=1, ..., n)
(3.20b)
lim
N →∞
bi and summing over i = 1,K ,n :
∑ i =1
bi2
(j=2, ..., m)
(3.20c)
( t k + i )Z d ( t k + i )T ≥ λ I
(I is the identity matrix)
i =0
The c dj ’s are in turn obtained simply multiplying (3.18b) by
d ij
d
2) The state vector Φ d ( t ) generated by (3.17) (or equally (3.10)) possesses the weak persistent excitation property i.e. there exists a real β>0, such that:
W1d ’s). Then, using (3.18b) one gets the bi ’s:
i
(3.22a)
]
That is, the coefficient c1d is perfectly known (because so is
i =1 n
)
v~ d ( t − 1 )Pmd ( v~ d ( t − 1 )) ...v~ d ( t − n )Pmd ( v~ d ( t − n )) (3.22b)
This question is answered making use of available properties. First, it readily follows from (3.11c) and (3.14b):
c =
(3.21)
v~ d ( t − 1 )P1d ( v~ d ( t − 1 )) ... v~ d ( t − n )P1d ( v~ d ( t − n )) ...
we obtain estimates of ( bi , c dj ) from those of µijd ?
∑b µ
( i )~ y ( i )T
[
As Θ d comes in linearly, equation (3.19c) turns out to be an adequate parameterization to get estimates of the parameters ai and µijd , using a least squares type algorithm. But, how can
d j
i =1
d
Z d ( t )T = − z d ( t − 1 )... − z d ( t − n )
v~ d (t − 1) Pmd (v~ d (t − 1))...v~ d (t − n) Pmd (v~ d (t − n))
n
N
∑Φ
Introduce the undisturbed regression vector:
(3.19a)
v~ d (t − 1) P1d (v~ d (t − 1)) ... v~ d (t − n) P1d (v~ d (t − n)) ...
bi = µ id1 / c1d = − µ id1 W1d
1 ( i )Φ d ( i )T N
(
Φ d (t ) T = [ − ~y (t − 1)... − ~y (t − n)
[
i =1
−1
d
A( q −1 ) z d ( t ) = B( q −1 )P d v~ d ( t )
Equation (3.18a) is given the linear regression form:
Θ d = a1 ... a n
N
∑Φ
This estimator will now be shown to be consistent. Let us consider the undisturbed version of model (3.17):
(3.18b)
~ y (t) = Φ d (t)T Θ d + ξ ( t )
(3.20d)
3.3 Parameter estimation
v~ d (t − i ) Pjd v~ d (t − i ) + ξ (t ) (3.18a)
µ ijd = bi c dj (i = 1,..., n; j = 1,..., m )
A(1) d y0 B (1)
1 N
N
∑Φ
d
( t )Φ d ( t )T > β I
(w.p.1)
i =1
3) The estimator (3.21) is consistent i.e. Θˆ d ( N ) → Θ d as N → ∞ (w.p. 1) Proof. The proof is omitted due to the limitation on the paper length. Its main ingredients can be found in (Chaoui et al., 2006).
346
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009 tk = k T
( k = 0 ,1,2K );
τ1 = 2n ;
τ j = ( 2 + j)n
Estimation of the bi ’s and the c j ’s. Relations (3.20b-c)
where
suggest the following estimators:
( j = 2 , ..., m ). The resulting limit cycle, denoted Alc , is shown in (Fig. 5). The rest of identification scheme development follows closely Section 3. The detailed design and analysis is omitted due to space limitation.
bˆi ( N ) = − µˆ i 1 ( N )W1d
(i=1, ..., n)
(3.24a)
(j=2, ..., m)
(3.24b)
n
∑ bˆ ( N )µˆ ( N ) i
cˆ dj ( N ) =
u
ij
i =1
n
∑ bˆ ( N ) 2 i
Cd
uM
i =1
(V
Proposition 3.2. The estimators (3.24a-b) are consistent. Proof. Comparing (3.24a) and (3.20b) yields ˆb ( N ) − b = (µ − µˆ ( N ))W d . The consistency of (3.24a) is i i i1 i1 1 then a direct consequence of the consistency of Θˆ d ( N )
(because this in particular ensures that µi 1 − µˆ i 1 ( N ) converges in probability to zero as N → ∞ ). The consistency of (3.24b) is proved similarly using the fact that both µˆ i 1 ( N ) and bˆi ( N ) are consistent estimators Estimation
of
the
points
(V
d j
)
Aˆ (q −1 , N ), Bˆ (q −1 , N )
A(q −1 ), B(q −1 )
are
the
estimates
by the previously obtained parameter estimates aˆ i ( N ), bˆi ( N ) . Comparing (3.25) and (3.20d) implies: Fˆ jd ( V jd , N ) − F d ( V jd ) =
(
)
hm
ˆ ( 1, N ) d A A( 1 ) d . y0 ( N ) − y0 B( 1 ) Bˆ( 1, N )
The consistency of (3.25) then follows from the consistency of cˆ dj ( N ), aˆ i ( N ), bˆi ( N ) and y 0d ( N ) . 4. IDENTIFICATION OF THE LINEAR SUBSYSTEM AND SWITCHED-HYSTERESIS ASCENDANT BORDER An identification scheme is developed that determines the linear subsystem parameters ( ai , bi ) as well as m points on the hysteresis ascendant border. There is a strong symmetry between the structure of this identification scheme and the one described in section 3. A key ingredient of the new scheme is an ascendant periodic experiment (APE) involving the periodic input signal v a ( t ) , with period T = n(m + 3 ) , defined as follows: V a for t = tk + τ j ; j = 1, ..., m v a(t) = j hm otherwise
)
a V2a V1 = hvM
V ma
Fig. 5. Movement of the working point on the limit cycle Alc when using the periodic input signal (4.1) 5. CONCLUSION
REFERENCES
of
induced
W jd cˆ dj ( N ) − c dj +
,C a (V jd+1 )
)
Hammerstein systems identification is considered is presence of memory nonlinearities of the switched-hysteresis type. The nonlinearity borders may be discontinuous, noninvertible and of unknown structure.
Aˆ (1, N ) d Fˆ d (V jd , N ) = cˆ dj ( N )W jd + y 0 ( N ) (3.25) Bˆ (1, N )
where
d j +1
,C a(V ja )
Ca
um
, F d (V jd ) , j = 2,..., m .
Equation (3.20d) suggests for F d ( V jd ) the estimator:
(V
Satellite
a j
Bai, E.W. (2002).‘Identification of linear systems with hard input nonlinearities of known structure’. Automatica, vol 38, pp. 853-860. Cerone, V., Regruto, D. (2007) ‘Bounding the parameters of linear systems with input backlash’. IEEE Trans. Automatic Control, vol 52(3), pp.531–536. Chaoui, F.Z., Giri, F., Rochdi, Y., Haloua, M., Naitali, A.(2005). ‘System identification based on Hammerstein model’. International Journal of Control, vol. 78, N° 6, pp. 430-442. Giri, F., Chaoui, F.Z., Rochdi, Y. (2002). ‘Interval excitation through impulse sequences. A technical lemma’. Automatica, vol 38, pp. 457-465. Giri, F., Rochdi, Y., Chaoui, F.Z., Brouri. A. (2008a). ‘Identification of Hammerstein systems in presence of hysteresis-backlash and hysteresis-relay nonlinearities’. Automatica, vol. 44, pp 767-775. Giri, F., Rochdi, Y., Elayan, E., Brouri, A., Chaoui, F.Z. (2008b) ‘Hammerstein systems identification in presence of hysteresis-backlash nonlinearity’. IFAC World Congress, Seoul, pp. 7859-7864. Greblicki, W., M. Pawlak (2008). ‘Nonparametric System Identification’. Cambridge University Press. 2008. ISBN-10: 0521868041.
(4.1)
347