Progressively better estimates of the domain of attraction for nonlinear systems via piecewise Takagi–Sugeno models: Stability and stabilization issues

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.1 (1-23)

Available online at www.sciencedirect.com

ScienceDirect Fuzzy Sets and Systems ••• (••••) •••–••• www.elsevier.com/locate/fss

Progressively better estimates of the domain of attraction for nonlinear systems via piecewise Takagi–Sugeno models: Stability and stabilization issues ✩ Temoatzin González, Miguel Bernal ∗ Sonora Institute of Technology, 5 de febrero 818 Sur, Col. Centro, C.P. 85000, Cd. Obregón, Mexico Received 23 September 2014; received in revised form 11 November 2015; accepted 13 November 2015

Abstract This paper introduces a novel approach for stability analysis and controller design of nonlinear models, both continuous- and discrete-time, based on an exact piecewise Takagi–Sugeno representation which generalizes the sector nonlinearity methodology. The motivation behind this work lies on the fact that piecewise-induced Takagi–Sugeno representations have smaller variations than single ones covering the whole modeling area, thus increasing the chances of finding piecewise solutions for stability analysis or controller synthesis. Moreover, this relaxation is shown to produce progressively better estimations of the domain of attraction via easy-to-implement algorithms that demand much less computational effort than other recent approaches in the literature. Examples are provided to highlight the advantages of the contributions over the existing methods. © 2015 Elsevier B.V. All rights reserved.

1. Introduction Originally intended to relate rule-based fuzzy representations with traditional control design, Takagi–Sugeno (TS) models first appeared in [1]; they were defined as a nonlinear convex blending of linear models where convex functions emerged from linguistic rules and linear models arose from linearization [2]: they were therefore approximate. Later on, the sector nonlinearity methodology allowed a nonlinear model to be exactly rewritten as a TS one in a compact set C of the state space, by capturing the model nonlinearities in membership functions (MFs) which held the convex-sum property [3]: this work adopts this approach. Thus, when exactly modeling a nonlinear system (a) a TS model is no longer an approximation, but an exact algebraic rewriting of the original one, regardless of where the state lies, (b) the  MFs hi are no longer state-independent entities belonging to a standard simplex  = hi ∈ R : 0 ≤ hi , i hi = 1 , ✩ This work has been supported by the Project Ciencia Básica SEP-CONACYT 168406 and the Project PROMEP ITSON-CA-18 1P4270406009 103.5/13/7283 Fortalecimiento. * Corresponding author. E-mail address: [email protected] (M. Bernal).

http://dx.doi.org/10.1016/j.fss.2015.11.010 0165-0114/© 2015 Elsevier B.V. All rights reserved.

JID:FSS AID:6943 /FLA

2

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.2 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

but state-dependent functions whose usefulness (i.e., whose convexity) is strictly limited to the modeling area C. This consideration is important since most of the TS-based results for stability analysis and controller design disregard their inherent locality, perhaps due to the fact that their final conditions do not take into account the MFs and look therefore alike whether they come from an exact or an approximate approach, i.e., they are linear matrix inequalities (LMIs) [4,5] whose feasibility can be efficiently tested via commercially available software [6]. Stability and stabilization conditions in the original TS-LMI framework are only sufficient [7,8]. One source of conservativeness is the way the MFs are removed from convex nested sums in order to obtain LMI conditions, which is called the co-positivity problem and has been given several answers [9–12], all of which disregard any state-dependence of the MFs. Another source of pessimism is the sort of Lyapunov function involved, which instead of quadratic can be fuzzy (reproducing the convex structure of the TS model) [13–18], or piecewise (for TS models where not all their linear components are simultaneously activated, thus inducing state-space partitions) [19–22]: the former is obliged to deal with the time-derivative of the MFs for continuous-time models, the latter is inapplicable for exact TS models where all the MFs are simultaneously activated. The use of more general classes of convex models has been also explored as a way to relax the aforementioned conservativeness: descriptor models [23,24] face the same problems of Lyapunov function choice; polynomial ones [25,26] may lead to better results at a great cost in theoretical complexity, numerical burden, and limited applicability. This paper considers that there is a compromise between complexity and quality of results that may have escaped from the picture above when stability and stabilization of nonlinear systems is investigated via TS models. To see this, consider a nonlinear system whose exact TS model is valid in C. Results based on a piecewise Lyapunov function Vp (x) (PWLF) associated to this system are valid for the largest region D = {x : Vp (x) ≤ γ } in C, while those based on a fuzzy Lyapunov function Vn (x) (FLF, also known as non-quadratic) are valid for the largest region  D= {x : Vn (x) ≤ γ } in C as long as a third set of bounds on the time-derivatives of the MFs hi hold: S = {x : h˙ i (x) ≤ φi } [27]. Therefore, piecewise approaches may lead to larger feasible regions. With this reasoning in mind, this paper has the following objectives: 1. Modeling: Since an exact representation should be available for analyzing nonlinear systems via TS models, this implies the use of the sector nonlinearity approach, but the TS models thus produced have all its linear consequents simultaneously activated, which precludes them from being treated under the piecewise approach. Efforts to overcome this problem have been made in [28,29] and [30] by mixing PWLFs with FLFs: they have produced very involved conditions where the continuity problem has not been addressed, except in [30] where the continuity is treated as in [19]. In [31] a generalization of the sector nonlinearity approach has been provided; it produces an exact piecewise TS (PWTS) representation of a nonlinear model: here, we pursue this investigation. 2. Stability analysis: Extend the methodology in [19] to exact PWTS models for stability analysis of both continuous- and discrete-time nonlinear systems showing (a) how an increasing number of partitions leads to progressively better feasibility sets which advantageously compare – both numerically and qualitatively – with other recent approaches such as [32–34] and (b) how the domain of attraction (DA) of an equilibrium point of a nonlinear system can be estimated with an arbitrary degree of accuracy via PWTS representations and an adequate algorithm, also performing better over recent similar approaches such as [26,27,35]. 3. Stabilization (controller design): Extend the methodologies in [21,36] to (a) determine a piecewise control law u(t) for nonlinear models via PWTS representations such that the closed-loop equilibrium point is locally asymptotically stable, both for continuous- and discrete-time systems, and (b) get progressively better invariant subsets of the DA that fit into the modeling region C improving over recent results concerned with the closed-loop DA [27,37]. Take into account that belonging to the modeling region C is an imposed condition in this sort of controller design, since the system thus modeled has physical constraints precisely reflected in the states belonging to C, whereas for stability of an autonomous system the DA may easily go beyond the modeling area. This article is organized in four more sections. Section 2 provides some notation, properties, and a generalization of the sector nonlinearity approach that allows a PWTS representation to be obtained from a nonlinear model on a compact set of its state space. Based on this representation, stability analysis is performed in Section 3 for continuousas well as discrete-time systems via PWTS models to the extent of providing progressively better estimations of the DA of the equilibrium point. Section 4 considers the stabilization of nonlinear systems via PWTS representations; the algorithm for estimation of the DA is adapted to controller design. Illustrative examples and comparisons are

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.3 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

3

provided along the paper to stress the advantages and simplicity of the proposed methodology over a number of techniques available in the literature. Section 5 gathers some concluding remarks. 2. A piecewise generalization of the sector nonlinearity approach Throughout this work, the following notation will be adopted: 

A Y

  (∗) A = B Y

YT B

 (1)

A + (∗) = A + AT .

(2)

Should a matrix expression be involved with symbols “>” and “<”, they will stand for positive and negativedefiniteness, respectively; “” and “≺” will stand for element-wise positive and negative, respectively. Arguments will be omitted when convenient. Consider the following affine-in-control nonlinear model: δx(t) = f (x(t))x(t) + g(x(t))u

(3)

where δx(t) = x(t) ˙ for continuous-time systems and δx(t) = x(t + 1) for discrete-time systems, x(t) ∈ Rn being the state vector, u(t) ∈ Rm being the input vector, f (·) : Rn → Rn×n and g(·) : Rn → Rn×m being smooth nonlinear vector fields, which are assumed to be bounded in a compact set C of the state space including the origin. Let Ck ⊂ C, k ∈ {1, 2, . . . , q} be a collection of compact subsets partitioning the modeling compact C of nonlinear model (3). The key idea behind the proposed approach is to apply the sector nonlinearity technique to each compact subset Ck , as to obtain q different TS representations of the nonlinear model (3). Each of these TS representations will preserve all the information in (3) since they all rewrite the original system, but the polytope associated with each convex representation may vary since it depends on the extreme values of the non-constant terms of f (x(t)) and g(x(t)) in Ck ; these variations give room to improvement as will be proven in the next section. Using the sector nonlinearity approach to represent nonlinear model (3) in q compacts yet preserving its unity, requires a proper notation. To this end, consider the non-constant terms of f (x(t)) and g(x(t)) in Ck as zjk (x(t)) ∈  zkj , zkj ; these vectors correspond to the premise or scheduling vector in traditional TS modeling. The fact that a TS representation of a nonlinear model is not unique comes from this stage, since the premise vector can be chosen in different ways [4]. The following weights can be constructed taking into account these bounds as jk

w0 (·) =

zkj − zjk (·) zkj − zkj

,

jk

jk

w1 (·) = 1 − w0 (·).

(4)

Consequently, r membership functions (MFs) for each of the q compact subsets Ck can be defined as hki (z(t)) =

p

j =1

jk wij zj ,

(5)

p , i ∈ {0, 1}, k ∈ {1, 2, . . . , q}. For each compact with i ∈ {1, 2, . . . , r}, i = ip × 2p−1 + . . . + i2 × 2 + i1 + 1, r = 2 j   subset Ck , the MFs thus defined hold the convex sum property, i.e., ri=1 hki · = 1, 0 ≤ hki · ≤ 1, k ∈ {1, 2, . . . , q}. Notice that the number of non-constant terms captured by the MFs hki does not vary from one compact subset to another. Taking into account that convex sums of matrix expressions ϒi , i ∈ {1, 2, . . . , r}, are usually written shortly as r  ϒz = hi (z(t)) ϒi , nonlinear model (3) can be written as the following PWTS model [31]: i=1

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.4 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

4

 δx t = Akz x (t) + Bzk u (t) , x(t) ∈ Ck , k ∈ {1, 2, . . . , q} ⎧r   h1i (z(t)) A1i x(t) + Bi1 u(t) = A1z x(t) + Bz1 u(t), x t ∈ C1 ⎪ i=1 ⎪ ⎪r   2 ⎪ 2 2 2 2 ⎨ x t ∈ C2 i=1 hi (z(t)) Ai x(t) + Bi u(t) = Az x(t) + Bz u(t), = . (6) ⎪.. ⎪ ⎪ ⎪   q ⎩r q q q q i=1 hi (z(t)) Ai x(t) + Bi u(t) = Az x(t) + Bz u(t), x t ∈ Cq   with Aki = f (x(t))hk =1 , Bik = g(x(t))hk =1 , i ∈ {1, 2, . . . , r}, k ∈ {1, 2, . . . , q}. i i It is important to stress that the original nonlinear model (3) is equivalent to the PWTS model (6), i.e., the PWTS model is not an approximation. On the other hand, note that this does not imply that the matrices in each submodel are the same, i.e., in general Aki = Ali , ∀k = l. As stated before, this methodology generalizes the sector nonlinearity approach in [3], since the latter corresponds to (6) with q = 1 (i.e., no partition at all, but one single modeling area C = C1 ). Example 1. To illustrate the aforementioned generalization of the sector nonlinearity approach, consider the following nonlinear model     1 −3.5 − 1.5a sin x1 (t) −4 x(t) + u(t) (7) x(t) ˙ = x2 (t) 9.5 − 10.5b sin x1 (t) −2 where the compact set C = {x : |x1 (t)| ≤ π, |x2 (t)| ≤ 2, } is the modeling area, and a, b ∈ Rn are known constants. Nonlinearities in (7) are z1 = sin x1 and z2 = x2 (p = 2), which are bounded as z1 = −1 ≤ z1 ≤ z1 = 1 and z2 = −2 ≤ z2 ≤ z2 = 2 in C. Therefore, this model can be rewritten in the TS form through the ordinary sector nonlinearity approach which corresponds to the methodology described above with q = 1. The corresponding TS model is: x(t) ˙ = Az x(t) + Bz u(t) =

r  i=1

hi (z(x(t))) (Ai x(t) + Bi u(t))

(8)

       −3.5 − 1.5az1 −4 −3.5 − 1.5az1 −4 1 1 , A 3 = A4 = , B1 = B3 = , B2 = B4 = , z2 z2 9.5 − 10.5bz1 −2 9.5 − 10.5bz1 −2 z1 − sin x1 2 z2 − x2 1 w01 = , w0 = , w = 1 − w01 , w12 = 1 − w02 , h1 = w01 w02 , h2 = w01 w12 , h3 = w11 w02 , h4 = w11 w12 . z1 − z1 z 2 − z2 1 Now consider a partition of C in two subsets (q = 2): C1 = {x : −π ≤ x1 (t) ≤ 0, |x2 (t)| ≤ 2} and C2 = {x : 0 ≤ x1 (t) ≤ π, |x2 (t)| ≤ 2}. Notice that the partition occurs only in x1 ; each subset contains the origin. Obviously, nonlinearities in (7) remain the same, but the bounds of those depending on x1 in C1 and C2 may have changed, i.e., z11 = −1 ≤ z11 ≤ z11 = 0 in C1 and z21 = 0 ≤ z12 ≤ z21 = 1 in C2 . Hence, applying the sector nonlinearity technique in each subset, the following PWTS model is obtained:      A1z x(t) + Bz1 u(t) = ri=1 h1i (z(x(t))) A1i x(t) + Bi1 u(t) , x t ∈ C1   x˙ t = (9)  A2z x(t) + Bz2 u(t) = ri=1 h2i (z(x(t))) A2i x(t) + Bi2 u(t) , x t ∈ C2         −3.5 + 1.5a −4 −3.5 −4 1 1 1 1 1 1 1 1 1 1 where A1 = A2 = , A3 = A4 = , B1 = B3 = , B2 = B4 = , A21 = 9.5 + 10.5b −2 9.5 −2 −2 2         −3.5 −4 −3.5 − 1.5a −4 1 1 , A23 = A24 = , B12 = B32 = , B22 = B42 = , w011 = − sin x1 , A22 = 9.5 −2 9.5 − 10.5b −2 −2 2 2 − x2 2 − x2 w021 = , w111 = 1 − w011 , w121 = 1 − w021 , w012 = 1 − sin x1 , w022 = , w112 = 1 − w012 , w122 = 1 − w022 , 4 4 h11 = w011 w021 , h12 = w011 w121 , h13 = w111 w021 , h14 = w111 w121 , h21 = w012 w022 , h22 = w012 w122 , h23 = w112 w022 , h24 = w112 w122 . T  Notice that, since the partition has been made only in x1 and g(x) = 1 x2 depends only on x2 , the set of Bik , i ∈ {1, . . . , 4} remains unchanged for k = 1 and k = 2. This model can be rewritten in a different PWTS form by dividing the compact set C into four subsets (q = 4): C1 = {x : −π ≤ x1 (t) ≤ 0, −2 ≤ x2 (t) ≤ 0}, C2 = {x : −π ≤ x1 (t) ≤ 0, 0 ≤ x2 (t) ≤ 2}, C3 = {x : 0 ≤ x1 (t) ≤ π, −2 ≤ x2 (t) ≤ 0}, and C4 = {x : 0 ≤ x1 (t) ≤ π, 0 ≤ x2 (t) ≤ 2}. Therefore, in each partition, the local TS representation is given by: 

where A1 = A2 =

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.5 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

x(t) ˙ =

r 

  hki (z(t)) Aki x(t) + Bik u(t) ,

5

x(t) ∈ Ck

(10)

i=1

where zk1 = max(sin x1 ), zk1 = min (sin x1 ), zk2 = max(x2 ), and zk2 = min (x2 ), w01k =

x∈Ck w11k = 1 − w01k ,  −3.5 − 1.5azk1 9.5 − 10.5bzk1

, w02k =

zk2 − x2

,

x∈Ck x∈Ck zk1 − zk1 zk2 − zk2 1k 2k k 1k 2k k 1k 2k k 1k 2k k = w0 w0 , h2 = w0 w1 , h3 = w1 w0 , h4 = w1 w1 , A1 = Ak2 =      −3.5 − 1.5azk1 −4 k = B k = 1 , B k = B k = 1 , k ∈ {1, . . . , 4}. , B 1 3 2 4 zk2 zk2 9.5 − 10.5bzk1 −2

x∈Ck − w02k , hk1 

w12k = 1  −4 , Ak3 = Ak4 = −2

zk1 − sin x1

3. Stability analysis based on exact PWTS models 3.1. The continuous-time case   Stability analysis of PWTS models naturally call for PWLFs of the form V x(t) = Vk x(t) for x(t) ∈ Ck , is giving conditions that translate Vk x(t) > 0 for x(t) = 0, k ∈ {1, 2, . . . , q}. An important issue  concerning  PWLFs into LMIs to guarantee its continuity, i.e., they must satisfy Vk x(t) = Vl x(t) , ∀x(t) ∈ (Ck ∩ Cl ). While some works simply assume continuity [28,29], others enforce it via polyhedral partitions [19]. The next developments adopt the latter approach for PWTS models. Notation in [19] will be adapted to the current approach. To this end, let K0 be the set of indexes of those compact subsets Ck that include the origin; otherwise, the indexes belong to the set K1 . The definitions    k  x(t) Ai 0 , x(t) ¯ = , A¯ ki = 1 0 0 i ∈ {1, 2, . . . , r}, k ∈ K1 , allow the PWTS model in (6) to be rewritten as: x(t) ˙ = Akz x(t), ˙¯ = A¯ kz x(t), ¯ x(t)

x(t) ∈ Ck , x(t) ∈ Ck ,

k ∈ K0 k ∈ K1 .





(11)

Polyhedral partitions of the state-space have the form E¯ k x¯ 0, x ∈ Ck , k ∈ {1, . . . , q} with E¯ = Ek ek such that ek = 0 for k ∈ K0 . Then, the following result directly arises from applying the piecewise stability analysis in [19] to the PWTS representation above: Theorem 1. If there exist symmetric matrices T , Uk 0, and Wki 0 such that  Pk = FkT T Fk , Pk − EkT Uk Ek > 0,  k T (12) k ∈ K0 , Ai Pk + Pk Aki + EkT Wki Ek < 0,  P¯k = F¯kT T F¯k , P¯k − E¯ kT Uk E¯ k > 0,  k T k ∈ K1 , (13) A¯ i P¯k + P¯k A¯ ki + E¯ kT Wki E¯ k < 0,   for i ∈ {1, 2, . . . , r}, where F¯ = Fk fk with fk = 0 for k ∈ K0 , satisfying F¯k x¯ = F¯l x¯ for x ∈ (Ck ∩ Cl ), k, l ∈ {1, 2, . . . , q}; then x t tends to zero exponentially for every continuous differentiable piecewise trajectory in C = q k=1 Ck satisfying (11). The corresponding piecewise Lyapunov function (PWLF) is given by: ⎧ T x ∈ Ck , k ∈ K0 ,  ⎨  xT Pk x,   (14) V x = x x ⎩ , x ∈ Ck , k ∈ K1 . P¯k 1 1 Proof. It follows directly from Theorem 1 in [19], provided that r matrices Aki or A¯ ki are always defined for each region Ck , thus explaining the inclusion of index set i ∈ {1, . . . , r} in (12) and (13). 2 Remark 1. It is worth noticing that stability conclusions derived from Theorem 1 are applicable to the original nonlinear model (3) because the PWTS representation (6) is algebraically equivalent to it; nevertheless, the convex

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.6 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

6

Fig. 1. (a) Comparison: “2” for quadratic case, “+” for Theorem 1 with 2 divisions, “×” for Theorem 1 with 4 divisions. (b) Comparison Example 1 [30]: “+” for Theorem 3 from [16], “◦” for Theorem 2 from [30], “×” for Theorem 1 with 16 divisions.

sum property only holds within the compact set C, because some MFs hi (z(t)) become negative or greater than one outside of C. This means that results in [19] are no longer limited to TS models with operating and interpolation regimes which can only approximate nonlinear systems. Remark 2. Matrices E¯ k , F¯k , k ∈ {1, 2, . . . , q} are not unique. A systematic procedure for their construction is described in [19] and [38]; additionally, the computational toolbox PWLTOOL can automatically generate them from the partition polyhedral boundaries [39]. Example 2. (a) Consider the nonlinear model in Example 1 with u(t) = 0, a ∈ [−4, 4], and b ∈ [−1.5, 1.5]. If the single TS model (8) is used for ordinary quadratic stability test, the feasible set of solutions marked with (2) in Fig. 1(a) is obtained. When the 2-region PWTS model (9) is used along with the piecewise stability analysis in Theorem 1, a larger feasible set of solutions is obtained as can be appreciated in Fig. 1(a) with marks (+). Matrices Fk , Ek , k ∈ {1, 2}, guaranteeing the continuity of the PWLF and allowing to implement the S-procedure relaxation, respectively, were obtained from the PWLTOOL [39]: ⎡ ⎤ ⎡ ⎤  T  T 0 0 1 0 −1 1 F1 = ⎣ 1 0 ⎦ , F 2 = ⎣ 1 0 ⎦ , E 1 = , E2 = . 0 0 0 1 0 1 Consider now that the modeling area C is further split into the following compact sets: C1 = {x : −π ≤ x1 (t) ≤ −π/4, |x2 (t)| ≤ 2} , C2 = {x : −π/4 ≤ x1 (t) ≤ 0, |x2 (t)| ≤ 2} , C3 = {x : 0 ≤ x1 (t) ≤ π/4, |x2 (t)| ≤ 2} , C4 = {x : π/4 ≤ x1 (t) ≤ π, |x2 (t)| ≤ 2} . Hence, a 4-region PWTS model arises: x(t) ˙ =

r 

hki (z(t))Aki x(t),

x(t) ∈ Ck ,

(15)

i=1

zk1 − sin x1



−3.5 − 1.5azk1 = 9.5 − 10.5bzk1

 −4 , Ak2 = −2

where = max(sin x1 ), = min (sin x1 ), = , =1− x∈Ck x∈Ck zk1 − zk1   −3.5 − 1.5azk1 −4 , k ∈ {1, . . . , 4}. As expected, this finer division of the state space produces a larger feasi9.5 − 10.5bzk1 −2 ble set of solutions for Theorem 1 than the aforementioned cases; this set is distinguished with (×) in Fig. 1(a). The following matrices Fk , Ek , k ∈ {1, . . . , 4} were also obtained using the PWLTOOL [39]: zk1

zk1

hk1

hk2

hk1 ,

Ak1

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.7 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

⎡

1 ⎢0 ⎢ F1 = ⎢ ⎢0 ⎣1 0  1 E1 = −1

⎤ 0 0.7854 0 0 ⎥ ⎥ 0 0 ⎥ ⎥, 0 0 ⎦ 1 0

⎡

⎡ ⎤ ⎤ 0 0 0 0 0 ⎢0 0⎥ ⎢ 1 0 −0.7854 ⎥ ⎢ ⎢ ⎥ ⎥ ⎢ ⎥ ⎥, 0 F3 = ⎢ 0 0 ⎥ , F4 = ⎢ ⎢0 0 ⎥ ⎣1 0⎦ ⎣1 0 ⎦ 0 0 1 0 1 0       −1 0 1 0 1 0 −0.7854 E2 = , E3 = , E4 = . −1 0 1 0 −1 0 1.5708

0 ⎢0 ⎢ F2 = ⎢ ⎢1 ⎣1 0

 0 1.5708 , 0 −0.7854

7

⎤ 0 0⎥ ⎥ 0⎥ ⎥, 0⎦ 1

⎡

(b) Example 1 in [16] and [30] has been considered for comparison with the proposed approach. To do so, the rule-based model has been rewritten as its nonlinear equivalent   x1 (sin x1 + 0.5 sin x2 − 3.5) − 4x2 x˙ = x1 (0.5(b − 1) + (1 + b)(0.2 sin x1 + 0.3 sin x2 )) − x2 ((2 + a)(0.3 sin x1 + 0.2 sin x2 ) + 0.5(2 − a)) and then remodeled with the piecewise methodology. As shown in Fig. 1(b), a 16-partition PWTS model along with conditions in Theorem 1 are enough to overcome the feasibility sets of previous approaches, even if they mix piecewise and non-quadratic methodologies as in [30]. 3.2. The discrete-time case Continuity of the PWLF in the discrete-time case is not required [40,20,41]; instead, the difficulty of the discretetime counterpart lies on the fact that the system trajectory can “jump” from one partition to another, not necessarily contiguous; therefore, all the combinations must be taken into account. To this end, some notation is introduced. Consider an autonomous discrete-time PWTS model: x(t + 1) = Akz x (t) ,

x(t) ∈ Ck ,

k ∈ {1, 2, . . . , q}

(16)

where q is the number of partitions. As before, the definitions:    k  x(t) Ai 0 x(t) ¯ = , A¯ ki = , k ∈ {1, 2, . . . , q} , 1 0 1 will be needed along with the following PWLF candidate:  T x Pk x, x ∈ Ck k ∈ K0 , V (x) = ¯ x ∈ Ck , k ∈ K1 , x¯ T P¯k x,

(17)

with K0 , K1 being sets of indexes grouping those regions which include the origin and those which do not, respectively; matrices Pk ∈ Rn×n : Pk = PkT > 0, k ∈ K0 and P¯k ∈ R(n+1)×(n+1) : P¯k = P¯kT > 0, k ∈ K1 . Additionally, consider the definitions   P 0 (18) , k ∈ K0 , Pˆk = PˆkT = k 0 0 which will come at hand to obtain LMI conditions guaranteeing V = V (x(t + 1)) − V (x(t)) < 0. As for the continuous-time case, the benefits of having locally-valid Lyapunov functions will be empowered by  the S-procedure. To that end, matrices of the form E¯ k x¯ 0, x(t) ∈ Ck , k ∈ {1, 2, . . . , q} with E¯ = Ek ek such that ek = 0 for k ∈ K0 , are defined. Theorem 2. If there exist matrices Pk = PkT , k ∈ K0 , P¯k = P¯kT , k ∈ K1 , Uk 0, and Wki 0 such that Pk − EkT Uk Ek > 0, k ∈ K0 , T P¯k − E¯ k Uk E¯ k > 0, k ∈ K1 ,  T  T Aki Pl Aki − Pk + EkT Wki Ek + Aki ElT Wli El Aki < 0,  T  T A¯ ki P¯l A¯ ki − P¯k + E¯ kT Wki E¯ k + A¯ ki E¯ lT Wli E¯ l A¯ ki < 0,

(19) (20) k, l ∈ K0 ,

(21)

k, l ∈ K1 ,

(22)

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.8 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

8

 T Pˆl A¯ ki − P¯k + E¯ kT Wki E¯ k + A¯ ki  T  T A¯ ki P¯l A¯ ki − Pˆk + E¯ kT Wki E¯ k + A¯ ki



A¯ ki

T

E¯ lT Wli E¯ l A¯ ki < 0,

k ∈ K1 , l ∈ K0 ,

(23)

E¯ lT Wli E¯ l A¯ ki < 0,

k ∈ K0 , l ∈ K1 ,

(24)

hold for i ∈ {1, 2, . . . , r}, then x(t) tends to zero exponentially for every trajectory in C = The corresponding piecewise Lyapunov function (PWLF) is given by (17).

q

k=1 Ck

satisfying (16).

Proof. Positive-definiteness as well as radial-unboundedness of V (x(t)) is guaranteed by (19) and (20) [19]. The time-difference V = V (x(t + 1)) − V (x(t)) is rewritten differently whether the departure point x(t) ∈ Ck corresponds to k ∈ K0 or k ∈ K1 , and if the arrival point x(t + 1) ∈ Cl corresponds to l ∈ K0 or l ∈ K1 . Therefore, there are four combinations reflected in LMIs (21)–(24). If the system “jumps” between the same kind of regions (i.e., k, l ∈ K0 or k, l ∈ K1 ), V < 0 if !  T V = x T (t + 1)Pl x(t + 1) − x T (t)Pk x(t) = x T (t) Akz Pl Akz − Pk x(t) < 0, (25) !  T V = x¯ T (t + 1)P¯l x(t ¯ + 1) − x¯ T (t)P¯k x(t) ¯ = x¯ T (t) A¯ kz P¯l A¯ kz − P¯k x(t) ¯ < 0. (26) The constructions of the constraint matrices Ek and E¯ k satisfying Ek x(t) 0, ¯ 0, E¯ k x(t)

∀x(t) ∈ Ck ,

k ∈ K0 ,

∀x(t) ∈ Ck ,

k ∈ K1 ,

with matrices Wki with nonnegative entries Wki 0 implies that x(t)T EkT Wki Ek x(t) ≥ 0, ∀x(t) ∈ Ck , T

x (t

+ 1)ElT Wli El x(t

k ∈ K0 ,

+ 1) ≥ 0, ∀x(t + 1) ∈ Cl ,

x(t) ¯ T E¯ kT Wki E¯ k x(t) ¯ ≥ 0, ∀x(t) ¯ ∈ Ck ,

l ∈ K0 ,

k ∈ K1 ,

x¯ T (t + 1)E¯ lT Wli E¯ l x(t ¯ + 1) ≥ 0, ∀x(t ¯ + 1) ∈ Cl ,

l ∈ K1 ,

for i ∈ {1, 2, . . . , r}, where sufficient conditions to guarantee (28) and (30) are  T x T (t) Aki ElT Wli El Aki x(t) ≥ 0, ∀x(t) ∈ Ck , k, l ∈ K0 ,  T x¯ T (t) A¯ ki E¯ lT Wli E¯ l A¯ ki x(t) ¯ ≥ 0, ∀x(t) ¯ ∈ Ck , k, l ∈ K1 ,

(27) (28) (29) (30)

(31) (32)

for i ∈ {1, 2, . . . , r}, respectively. Therefore, LMIs (21) using S-procedure with (27) and (31) guarantee (25), and (22) using S-procedure with (29) and (32) guarantee (26). As for “jumps” between regions which do not belong to the same kind, two cases are distinguished: if k ∈ K1 , l ∈ K0 , then V < 0 if !  T V = x¯ T (t + 1)Pˆl x(t ¯ + 1) − x¯ T (t)P¯k x(t) ¯ = x¯ T (t) A¯ kz Pˆl A¯ kz − P¯k x(t) ¯ < 0. (33) As in the previous case, the following inequalities hold: x(t) ¯ T E¯ kT Wki E¯ k x(t) ¯ > 0, ∀x(t) ¯ ∈ Ck , k ∈ K1 ,  T x¯ T (t) A¯ ki E¯ lT Wli E¯ l A¯ ki x(t) ¯ > 0, ∀x(t) ¯ ∈ Ck , k, l ∈ K0 .

(34) (35)

Thus, LMIs (23) using S-procedure with (34) and (35) guarantee (33). Similarly, if k ∈ K0 , l ∈ K1 , then V < 0 if !  T V = x¯ T (t + 1)P¯l x(t ¯ + 1) − x¯ T (t)Pˆk x(t) ¯ = x¯ T (t) A¯ kz P¯l A¯ kz − Pˆk x(t) ¯ < 0. (36)

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.9 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

9

Likewise as before, the following inequalities hold: x(t) ¯ T E¯ kT Wki E¯ k x(t) ¯ > 0, ∀x(t) ¯ ∈ Ck , k ∈ K0 ,  T ¯ > 0, ∀x(t) ¯ ∈ Ck , k, l ∈ K1 . x¯ T (t) A¯ ki E¯ lT Wli E¯ l A¯ ki x(t)

(37) (38)

Therefore, LMIs (24) using S-procedure with (37) and (38) guarantee (36), thus concluding the proof. Example 3. For the sake of comparison, consider the following discrete-time TS model:    (0.45 − 0.5a)x1 + 0.5a + 0.45 x1 0.4(x1 − 1) x(t + 1) = 0.5b − 0.45 + (0.45 + 0.5b)x1 x2 −0.4(x1 + 1)

(39)

Rewriting this model in the TS form within the compact set C = {x : |x1 (t)| ≤ 1, x2 (t) ∈ R}, the following TS model is obtained: x(t + 1) = Az x(t) =

r 

(40)

hi (z(x(t)))Ai x(t)

i=1

   1 − x1 a −0.8 0.9 0 , A2 = , h1 = where A1 = , h2 = 1 − h1 . In this example, the quadratic framework 0 −0.9 −0.8 b 2 does not have solution when the parameters a and b are in the interval a ∈ [−1, 1] and b ∈ [−1, 1], though simulations suggest that the system is stable [32]; therefore, the proposed piecewise approach comes at hand. Consider a partition of C in two subsets (q = 2): C1 = {x : −1 ≤ x1 (t) ≤ 0, x2 (t) ∈ R} and C2 = {x : 0 ≤ x1 (t) ≤ 1, x2 (t) ∈ R}. Hence, the following 2-region PWTS model arises: 

x(t + 1) =

r 

hki (z(t))Aki x(t),

x(t) ∈ Ck ,

(41)

i=1

where zk1 = max(x1 ), zk1 = min (x1 ), hk1 = x∈Ck



x∈Ck

zk1 − x1 zk1 − zk1

, hk2 = 1 − hk1 , and

 (0.45 − 0.5a)zk1 + 0.5a + 0.45 0.4(zk1 − 1) , −0.4(zk1 + 1) 0.5b − 0.45 + (0.45 + 0.5b)zk1   (0.45 − 0.5a)zk1 + 0.5a + 0.45 0.4(zk1 − 1) k A2 = , −0.4(zk1 + 1) 0.5b − 0.45 + (0.45 + 0.5b)zk1

Ak1 =

for k ∈ {1, 2}. In contrast with the quadratic case, Theorem 2 is feasible for a number of instances (a, b) of the PWTS representation above, whose partitions are described by matrices     E1 = −1 0 , E2 = 1 0 . These instances are shown in Fig. 2(a) with (•). Following the same methodology, a 4-region PWTS model was obtained with the following partitions: C1 = {x : −1 ≤ x1 (t) ≤ −0.5, x2 (t) ∈ R} , C3 = {x : 0 ≤ x1 (t) ≤ 0.5, x2 (t) ∈ R} ,

C2 = {x : −0.5 ≤ x1 (t) ≤ 0, x2 (t) ∈ R} ,

C4 = {x : 0.5 ≤ x1 (t) ≤ 1, x2 (t) ∈ R} ,

which can be described by the following matrices:     1 0 1 −1 0 ¯ , E2 = , E1 = −1 0 −0.5 −1 0



 1 0 E3 = , 1 0

E¯ 4 =



1 0 −1 0

 −0.5 . 1

This partition is the result of further division of the previous one; it is therefore expected to produce at least the same feasible set as the previous one when Theorem 2 is applied. This fact can be verified in Fig. 2(a) where the feasible set is indicated by (2).

JID:FSS AID:6943 /FLA

10

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.10 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

Fig. 2. Comparison Example 3: “•” for Th. 2 with q = 2, “2” for Th. 2 with q = 4, “×” for Th. 2 with q = 6, “+” for Th. 2 with q = 92, “◦” for maximum feasible set in Fig. 7 of [32].

Another partition refining the original 2-region PWTS model can be the following: C1 = {x : −1 ≤ x1 (t) ≤ −2/3, x2 (t) ∈ R} ,

C2 = {x : −2/3 ≤ x1 (t) ≤ −1/3, x2 (t) ∈ R} ,

C3 = {x : −1/3 ≤ x1 (t) ≤ 0, x2 (t) ∈ R} ,

C4 = {x : 0 ≤ x1 (t) ≤ 1/3, x2 (t) ∈ R} ,

C5 = {x : 1/3 ≤ x1 (t) ≤ 2/3, x2 (t) ∈ R} ,

C6 = {x : 2/3 ≤ x1 (t) ≤ 1, x2 (t) ∈ R} ,

with the following descriptive matrices:       1 0 1 1 0 2/3 −1 0 , E¯ 2 = , E3 = , E¯ 1 = −1 0 −2/3 −1 0 −1/3 −1 0       1 0 1 0 −1/3 1 0 −2/3 , E¯ 5 = , E¯ 6 = . E4 = 1 0 −1 0 2/3 −1 0 1 Note that this 6-region partition is not a refinement of the 4-region one above; this means that the new feasibility set does not necessarily contain that of the 4-region partition (though for this example, it does). Marks (×) in Fig. 2(a) show the feasibility area of the 6-partition PWTS model when Theorem 2 is applied: it is compared with former partitions (q = 2, 4). For comparison purposes, a 92-region PWTS model is considered dividing the modeling area C into 92 compact sets. In Fig. 2(a) the feasible area of the 92-partition PWTS model applying Theorem 2 is indicate by (+), where this finer division of the state space produces a larger feasible set of solutions than the above cases. For the sake of comparison, this feasible set is shown in Fig. 2(b) along with the maximum feasibility set in Fig. 7 of [32], where the most complex k-former-samples approach has been used with k = 3 and denoted with (◦). As illustrated, the proposed methodology based on exact PWTS models can outperform complex non-quadratic analysis which require remarkably higher computational effort [33,32]. 3.3. Estimation of the domain of attraction The goal of this section is to obtain the largest possible estimate of the DA of a nonlinear system with an asymptotically stable origin, no matter if it has more than one equilibrium point. The following definition of DA will be adopted: Domain of attraction [42] is the set = {x : ∀t ≥ 0, limt→∞ φ(t, x) = 0} for a given nonlinear system (7) (u(t) = 0) with its origin x = 0 being asymptotically stable and φ(t, x) being its solution. Unless nonlinearities are naturally bounded for the whole state space Rn of a nonlinear system, any exact TS representation – whether simple or piecewise – will be only locally valid for a compact set C0 of the state space. This is due to the fact that MFs – which depend on the system nonlinearities – will not stand for the convex-sum property outside the modeling area C0 , thus invalidating any Lyapunov-based analysis beyond this compact. Thus, the outermost

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.11 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

11

Fig. 3. An illustration of a 1-hypercube in R2 .

Lyapunov set 0 = V (x) ≤ γ within C0 is the maximum valid estimation of the DA through TS-LMI/Lyapunov-based methods. This is the approach in [17,43,27], where nested non-quadratic FLFs have been used. This method has two drawbacks: a rapidly growing computational complexity and the inclusion of a third set besides the modeling compact C0 and  the outermost Lyapunov set 0 in order to take into account the time-derivatives of the MFs, i.e., S = {x : h˙ i (x) ≤ φi }. Following the same principles, several results on Lyapunov-based estimation of the DA via sum of squares have recently appeared [44,45,26]. SOS-based estimation of the DA faces even more involved numerical calculations since locality requires adaptations of the Positivstellensatz theorem in order to add multipliers to check local positiveness. Some recent methods have also appeared that without a Lyapunov function can provide estimations of the DA based on exact TS or polynomial representations of a nonlinear system [37,35]. The results that follow belong to the first class of solutions, i.e., they are Lyapunov-based, though they employ PWLFs and exact PWTS models. The intuition behind the proposed method is that progressively better estimates can be found if the state space is iteratively partitioned as to asymptotically fit the DA. This partition should be made according to the states involved in the system nonlinearities that are captured in the MFs of any exact PWTS representation since the “inner” TS model within each partition will exhibit variations leading to bigger chances of feasibility. The algorithm begins by finding any proven subset of the domain of attraction C0 – e.g., by quadratic methods – and then using it as the initial set for piecewise partitions surrounding the proven subset: those “added” partitions that lead to feasible conditions under PWTS stability analysis remain; those who don’t are dismissed. A refinement of the partition follows along with an enlargement of the tested zone,  and the procedure is repeated iteratively. Let S = xi1 , xi2 , . . . , xis be the set of s states that are present in the nonlinearities of a given nonlinear system (3) with u(t) = 0, and Is = {i1 , i2 , . . . , is } the set of indexes referring to those states. From the definitions in Section 2, it is clear that the premise variables z(x(t)) and therefore the corresponding MFs hi (z(t)) depend exclusively on   these states. Let H = x ∈ Rn : xi0 − α/2 ≤ xi ≤ xi0 + α/2, i ∈ Is be the s-hypercube of length α ∈ R+ and center in xi = xi0 , i ∈ Is . In order to clarify these concepts, consider the model (7) with u(t) = 0. Clearly, nonlinearities depend solely on x1 , i.e., S = {x1 } and I1= {1}; based on this, Fig. 3 illustrates the corresponding 1-hypercube of length 1 and center x1 = x10 = 0, i.e., H = x ∈ R2 : |x1 | ≤ 1 . Let C0 be the s-hypercube of maximum length and center at the origin on which an exact TS model of a nonlinear system is proven to be quadratically stable, i.e., C0 = {x ∈ Rn : |xi | ≤ α0 , i ∈ Is } with α0 = max |xi | such that ∃P = i  P T > 0, P Aj + ATj P < 0, j = 1, 2, . . . , r, for δx(t) = Az x(t), i hi = 1, hi ≥ 0 modeled within C0 .

JID:FSS AID:6943 /FLA

12

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.12 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

Fig. 4. Description of Algorithm 1.

Algorithm 1. Starting from the compact C0 defined above – the maximum s-hypercube on which the nonlinear system is proven to be quadratically stable via a TS model – with length α0 , set i = 1 and α1 = α0 and perform the following steps: j

1. Consider s-hypercube compacts Ci of length αi indexed by j ∈ Ii and fully surrounding Ci−1 by sharing a side j or a vertex with it, such that Ci ∩ Cik = ∅ for j = k. 2. Set Ci = Ci−1 . j 3. For each j ∈ Ii test whether Ci ∪ Ci under its corresponding PWTS representation renders Theorem 1 (if continj uous) or 2 (if discrete) feasible or not. If yes, redefine Ci = Ci ∪ Ci . Go to the next j until all the indexes in Ii have been run. j 4. If Ci includes all the regions Ci , j ∈ Ii (i.e., all of them altogether with Ci−1 are feasible in the PWTS sense), then set αi+1 as the length of Ci ; otherwise, set αi+1 as a half of the length of any Ci−1 . 5. The current estimation of the DA is given by the outermost Lyapunov set V (x) ≤ γ in Ci with V (x) as the PWLF (14) (continuous-time case) or (17) (discrete-time case). If this estimate is satisfactory or the LMIs are rendered computationally intractable, terminate the algorithm; otherwise, set i = i + 1 and go to step 1 for a finer estimation of the DA. Fig. 4 illustrates the algorithm above for a 2nd-order system with nonlinearities in both states: x1 and x2 . It begins by calculating the maximum s-hypercube C0 on which the system is feasible and goes on to test a bigger area with partitions of the same size: those which can be added rendering feasible conditions in Theorem 1 (if continuous) or Theorem 2 (if discrete), remain (shadowed areas); otherwise, they are dismissed. The outermost Lyapunov set in this area is an estimation of the DA. If a finer estimation is required, the procedure is repeated by adding new sets around

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.13 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

13

Fig. 5. Differences between Algorithm 1 and Algorithm 2.

C1 (half the size of the previous one) and testing again for feasibility. The procedure is repeated up to numerical tractability or the desired degree of precision. A slight modification of the previous algorithm may help to improve its feasibility chances, though at the price of a higher computational cost. It is based on the following property: Lemma 3. qIf a nonlinear system (7) with u(t) = 0 is proven to be stable via an exact PWTS model with a given partition i=1 Ci , along with conditions in Theorem 1 (if continuous) or Theorem 2 (if discrete), so it will be for any q q q j j finer partition i=1 j i=1 Ci , such that Ci = j i=1 Ci . q Proof. If partition i=1 Ci proved feasible, it implies the existence of a PWLF (14) (or (17) for the discrete case) with q q j j a different Vi (x(t)) for each partition i ∈ {1, 2, . . . , q}. For a finer partition i=1 j i=1 Ci it suffices to take Vi = Vi for all j ∈ {1, 2, . . . , qi } to guarantee the feasibility of the new conditions since they will reduce to the original ones. Thus, step 1 in the algorithm above can be replaced by the following: j

1. Consider s-hypercube compacts Ci of length αi indexed by j ∈ Ii and fully covering and surrounding Ci−1 by j sharing a side or a vertex with it, such that Ci ∩ Cik = ∅ for j = k. This modification – referred in the sequel as Algorithm 2 – enables finer divisions within inner already-proven partitions, thus increasing the feasibility chances of the whole test at every step by providing the Lyapunov function additional flexibility. Fig. 4 is reproduced in Fig. 5 along with finer divisions as to point out the differences. Example 4. Consider again the nonlinear model in Example 1 with u(t) = 0, a = 1, and b = 1. The objective is to estimate the largest possible DA of the origin using Algorithm 1. This task begins by finding the maximum s-hypercube C0 on which the system is quadratically stable, which results in C0 = {x ∈ R2 : |x1 | ≤ 0.997} (C41 in Fig. 6(a)). For illustration purposes, Algorithm 1 runs only for four iterations, but it can go further if better estimations are needed up to the capabilities of the LMI solver; the system is indeed globally asymptotically stable. The improvement over estimates in [27] as well as others found therein can be checked in Fig. 6(b). Four system trajectories are also included to corroborate the validity of the Lyapunov curves. Example 5. Consider the following nonlinear model from [17]:     −0.2363 + 0.0985 x12 + x22 −1  x1 x(t) ˙ = x2 x12 + 2 −0.7097 + 0.3427 x12 + x22

(42)

Phase-plane simulations show that the system has an unstable limit cycle, which implies the DA is bounded by it. In order to obtain the largest possible estimate of the DA, Algorithm 2 comes at hand. As before, this starts by obtaining

JID:FSS AID:6943 /FLA

14

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.14 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

Fig. 6. Estimates of the DA for Example 4: Algorithm 1 (outermost black closed curve), Th. 1 in [46] (green closed curve), Th. 1 in [27] (purple closed curve), Th. 2 in [17] (blue closed curve), Th. 1 in [47] (red closed curve). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Fig. 7. Estimates of the DA for Example 5: limit cycle (outermost green closed curve), Algorithm 2 (black closed curve), Th. 2 in [17] (blue closed curve). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

the maximum s-hypercube C0 on which the system is quadratically stable; it has been found that this set corresponds to C0 = {x ∈ R2 : |xi | ≤ 0.9791}, i ∈ {1, 2}, which was used as an initial estimate of the DA. Fig. 7 plots the limit cycle (outermost green closed curve) and compares it with a first estimation of the DA given by one iteration of Algorithm 2 (black closed curve in Fig. 7(a)); then, an estimation of the DA given by 3 iterations of Algorithm 2 (192 partitions) is shown in Fig. 7(b) as a black closed curve and compared with the estimation provided in [17] (blue closed curve). 4. Controller design based on exact PWTS models This section is concerned with controller design of nonlinear systems based on an exact PWTS representation (6). The use of PWLFs for controller design of continuous-time TS models was first considered in [40], but the LMI conditions thereby found turned out to be erroneous [41]. A BMI solution was offered in [21], though BMI solvers at the time were not available. In all these cases, the TS model was assumed not to have all their linear models activated at once; unfortunately, this assumption does not hold for TS models built by using the sector nonlinearity approach [5]. As for the discrete-time case, the works of [48,49,22] provide LMI conditions for stabilization, though again they are restricted to TS models which do not have all their linear models simultaneously activated. These obstacles have been removed through the new methodology as shown in this section.

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.15 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

15

4.1. The continuous-time case We begin by rewriting TS model (6) as: x(t) ˙ = Akz x(t) + Bzk u(t), ˙¯ = A¯ kz x(t) ¯ + B¯ zk u(t), x(t)

x(t) ∈ Ck , x(t) ∈ Ck ,

k ∈ K0 , k ∈ K1 ,

with Az , A¯ z , K0 , K1 , x¯ having the same meanings as in the previous section, and B¯ ik = {1, 2, . . . , r}, k ∈ K1 . Consider the following piecewise PDC control law: u(t) = Lkz x(t), ¯ u(t) = L¯ kz x(t),

x(t) ∈ Ck , x(t) ∈ Ck ,

k ∈ K0 , k ∈ K1

(43) 

Bik

T

T 0

, i∈

(44)

where Lki ∈ Rm×n , i ∈ {1, . . . , r}, k ∈ K0 , and L¯ ki ∈ Rm×(n+1) , i ∈ {1, . . . , r}, k ∈ K1 . Substituting (44) in (43) gives the closed-loop PWTS model:  x(t) ˙ = Akz + Bzk Lkz x(t), x(t) ∈ Ck , k ∈ K0 , ˙¯ = A¯ kz + B¯ zk L¯ kz x(t), ¯ x(t) ∈ Ck , k ∈ K1 . x(t)

(45)

Theorem 4. If there exist symmetric matrices T , Uk 0, and Wki 0, and matrices Lki and L¯ ki such that " Pk = FkT T Fk , Pk − EkT U Ek > 0, (46) k ∈ K0 , Pk Aki + Pk Bik Lkj + (∗) + EkT Wki Ek < 0, " P¯k = F¯kT T F¯k , P¯k − E¯ kT U E¯ k > 0, (47) k ∈ K1 , P¯k A¯ ki + P¯k B¯ ik L¯ kj + (∗) + E¯ kT Wki E¯ k < 0,     hold for i, j ∈ {1, 2, . . . , r}, where E¯ = Ek ek with ek = 0 for k ∈ K0 , satisfying E¯ k x¯ 0 for x ∈ Ck ; F¯ = Fk fk  with fk = 0 for k ∈ K0 , satisfying F¯k x¯ = F¯l x¯ for x ∈ (Ck ∩ Cl ), k, l ∈ {1, 2, . . . , q}; then x t tends to zero exponenq tially for every continuous differentiable piecewise trajectory in C = k=1 Ck satisfying (11). Proof. It follows directly from proof of Theorem 1. Example 6. Consider the following nonlinear model in [50]:      0.5a (sin x1 + 1) − sin x1 + 1 2.5 sin x1 − 7.5 x1 0.5b (sin x1 + 1) − 0.5 (sin x1 − 1) u(t). x(t) ˙ = + 0.5 (sin x1 + 3) sin x1 + 1 x2 1.5 − 0.5 sin x1 (48) Taking z(x) = sin x1 , this model can be rewritten in a TS form within the compact C = {x : |xi | ≤ 2π }, i ∈ {1, 2}, obtaining the following TS model: x(t) ˙ =

2  i=1



where A1 =

hi (z(t)) (Ai x(t) + Bi u(t))

  2 −10 a , A2 = 2 0 1

(49)

     −5 1 b , B1 = , B2 = , h1 (z(t)) = 0.5 − 0.5 sin x1 (t), h2 (z(t)) = 1 − 2 1 2

h1 (z(t)). Consider the case where a = −20 and b = 4; ordinary quadratic conditions are unfeasible for this case. Thus, the proposed approach comes at hand. To this end, the modeling compact C is divided in 4 regions C1 = {x : −2π ≤ x1 ≤ −π, |x2 | ≤ 2π}, C2 = {x : −π ≤ x1 ≤ 0, |x2 | ≤ 2π }, C3 = {x : 0 ≤ x1 ≤ π, |x2 | ≤ 2π }, C4 = {x : π ≤ x1 ≤ 2π, |x2 | ≤ 2π }, on which the following 4-region PWTS model is defined: x(t) ˙ =

r  i=1

  hki (z(t)) Aki x(t) + Bik u(t) ,

x(t) ∈ Ck ,

(50)

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.16 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

16

Fig. 8. Results in Example 6.

where zk = max(sin x1 ), zk = min (sin x1 ), hk1 =

zk − sin x1

, hk2 = 1 − hk1 , and zk − zk      0.5a zk + 1 − zk + 1 2.5zk − 7.5 0.5a zk + 1 − zk + 1 k Ak1 = = , A 2 1.5 − 0.5zk zk + 1 1.5 − 0.5zk       0.5b z + 1 − 0.5 z − 1 0.5b (z + 1) − 0.5 (z − 1) , B2k = , B1k = 0.5 (z + 3) 0.5 z + 3 x∈Ck

x∈Ck

 2.5zk − 7.5 , zk + 1

for k ∈ {1, . . . , 4}. The following matrices Fk , Ek , k ∈ {1, . . . , 4} are obtained using the PWLTOOL [39]:         1 0 2π −1 0 1 0 1 0 −π E1 = , E2 = , E3 = E4 = , −1 0 −π −1 0 1 0 −1 0 2π ⎡ ⎡ ⎡ ⎡ ⎤ ⎤ ⎤ ⎤ −1 0 −π 0 0 0 0 0 0 0 ⎢ 0 0 0 ⎥ ⎢0 0⎥ ⎢0 0⎥ ⎢ −1 0 π ⎥ ⎢ ⎢ ⎢ ⎢ ⎥ ⎥ ⎥ ⎥ ⎢ ⎢ ⎢ ⎥ ⎥ ⎥ ⎥ F1 = ⎢ ⎢ 1 0 0 ⎥ , F2 = ⎢ 1 0 ⎥ , F3 = ⎢ 0 0 ⎥ , F4 = ⎢ 0 0 0 ⎥ . ⎣ 1 0 0 ⎦ ⎣1 0⎦ ⎣1 0⎦ ⎣ 1 0 0⎦ 0 1 0 0 1 0 1 0 1 0 Due to the fact that conditions in Theorem 5 are expressed as BMIs, the PENBMI toolbox in MATLAB [51] is used to obtain the following feasible solution:     L11 = −1.6871 −15.4875 −2.2889 , L12 = 8.6413 −10.5703 −1.2649 ,     L21 = 6.2962 −17.0926 , L22 = 7.7821 −16.7110 ,     L31 = −2.8254 −9.7310 , L32 = −1.0415 −9.2444 ,     L41 = −0.0441 −2.1607 2.6473 , L42 = −1.1205 −2.3192 1.9380 ,     2.0216 −2.1552 2.6414 −1.1572 P2 = 1 × 107 , P3 = 1 × 107 , −2.1552 4.5811 −1.1572 4.5811 ⎡ ⎡ ⎤ ⎤ 2.0086 −1.8445 −0.4567 0.7823 −0.8012 0.9194 0.9762 ⎦ , P4 = 1 × 107 ⎣ −0.8012 4.5811 −1.1181 ⎦ . P1 = 1 × 107 ⎣ −1.8445 4.5811 −0.4567 0.9762 −2.7412 0.9194 −1.1181 12.5708 Fig. 8(a) shows with q solid lines the boundaries of each subcompact Ck , k ∈ {1, . . . , 4}; some Lyapunov curve levels of (14) within C = k=1 Ck are shown in dashed lines. For the sake of comparison, the boundaries of the biggest region of attraction founded in [34] is shown with a dotted line: it is overlapped with the biggest region of attraction found via Theorem 5. In Fig. 8(b) it is shown the time evolution of the control law for the trajectories in Fig. 8(a).

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.17 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

17

4.2. The discrete-time case In this case, the PWTS model from (6) is rewritten as: x(t + 1) = Akz x(t) + Bzk u(t), x(t ¯ + 1) = A¯ kz x(t) ¯ + B¯ zk u(t),

x(t) ∈ Ck , k ∈ K0 x(t) ∈ Ck , k ∈ K1 , T  where Az , A¯ z , K0 , K1 , x, ¯ and B¯ ik = Bik T 0 have the same meanings as in the previous subsection. As before, the piecewise PDC control law (44) is considered. Additionally consider the definition   Lˆ kz = Lkz 0 , k ∈ K0 . The corresponding closed-loop PWTS model obtained by substituting (44) in (51) is given by:  x(t + 1) = Akz + Bzk Lkz x(t), x(t) ∈ Ck , k ∈ K0 ¯ x(t) ∈ Ck , k ∈ K1 . x(t ¯ + 1) = A¯ kz + B¯ zk L¯ kz x(t),

(51)

(52)

Theorem 5. If there exist matrices Pk = PkT , k ∈ K0 , P¯k = P¯kT , k ∈ K1 , Uk 0, and Wki 0, and matrices Lki , L¯ ki , ¯ kl , Hlk , and Glk of proper dimension such that Gkl , G Pk − EkT Uk Ek > 0, k ∈ K0 T ¯ ¯ ¯ k ∈ K1 Pk − Ek Uk Ek > 0, # k k k Gkl Ai + Gkl Bi Lj + (∗) − Pk + EkT Wki Ek Hkl Aki + Hkl Bik Lkj − GTkl # ¯ kl B¯ k L¯ k + (∗) − P¯k + E¯ T Wki E¯ k ¯ kl A¯ k + G G k i i j ¯T H¯ kl A¯ ki + H¯ kl B¯ ik L¯ kj − G kl # k k k ¯ ¯ ¯ ¯ ¯ ¯ Gkl Ai + Gkl Bi Lj + (∗) − Pk + E¯ kT Wki E¯ k ¯T H¯ kl A¯ ki + H¯ kl B¯ ik L¯ kj − G kl # k k k ¯ kl A¯ + G ¯ kl B¯ Lˆ + (∗) − Pˆk + E¯ T Wki E¯ k G k i i j ¯T H¯ kl A¯ k + H¯ kl B¯ k Lˆ k − G i

i

j

kl

(53) (∗) −Hkl − HklT + Pl + ElT Wli El (∗) T ¯ ¯ −Hkl − Hkl + P¯l + E¯ lT Wli E¯ l (∗) −H¯ kl − H¯ klT + Pˆl + E¯ lT Wli E¯ l (∗) T ¯ ¯ −Hkl − Hkl + P¯l + E¯ lT Wli E¯ l

(54)

$ < 0, k, l ∈ K0 ,

(55)

< 0, k, l ∈ K1 ,

(56)

< 0, k ∈ K1 , l ∈ K0

(57)

< 0, k ∈ K0 , l ∈ K1

(58)

$ $ $

hold for i ∈ {1, 2, . . . , r}, and k, l ∈ {1, 2, . . . , q}, then x(t) tends to zero exponentially for every trajectory in C =  q k=1 Ck satisfying (16). The corresponding piecewise Lyapunov function (PWLF) is given by (17). Proof. Consider the PWLF (17). As in the stability case, the time-difference V = V (x(t + 1)) − V (x(t)) is rewritten differently whether the departure point x(t) ∈ Ck corresponds to k ∈ K0 or k ∈ K1 , and if the arrival point x(t + 1) ∈ Cl corresponds to l ∈ K0 or l ∈ K1 . Therefore, there are four combinations reflected in LMIs (21)–(24). If the system “jumps” between the same kind of regions (i.e., k, l ∈ K0 or k, l ∈ K1 ), V < 0 if  T    x(t) −Pk 0 x(t) V = < 0, (59) x(t + 1) 0 Pl x(t + 1)  T    −P¯k 0 x(t) ¯ x(t) ¯ < 0. (60) V = x(t ¯ + 1) x(t ¯ + 1) 0 P¯l From expressions (27), (28), (29) and (30) we have  T  T    x(t) Ek W k E k x(t) 0 ≥0 x(t + 1) x(t + 1) 0 ElT Wl El   T  T   x(t) ¯ 0 E¯ k Wk E¯ k x(t) ¯ ≥0 x(t ¯ + 1) x(t ¯ + 1) 0 E¯ lT Wl E¯ l

∀x(t), ∈ Ck , k ∈ K0 , ∀x(t + 1), ∈ Cl , l ∈ K0 ,

(61)

∀x(t), ¯ ∈ Ck , k ∈ K1 , ∀x(t ¯ + 1), ∈ Cl , l ∈ K1 .

(62)

JID:FSS AID:6943 /FLA

18

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.18 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

The closed-loop models from (52) can be rewritten as:   "  k  x(t) ∀x(t), ∈ Ck , k ∈ K0 , k k =0 Az + Bz Lz −I ∀x(t + 1), ∈ Cl , l ∈ K0 , x(t + 1)   "   k x(t) ¯ ∀x(t), ¯ ∈ Ck , k ∈ K1 , k k ¯ ¯ ¯ =0 Az + Bz Lz −I x(t ¯ + 1) ∀x(t ¯ + 1), ∈ Cl , l ∈ K1 . Combining inequalities (61)–(64), with the previous equations by Finsler’s lemma,1 yields:    "   0 G  k −Pk + Ek Wk Ek < 0 k, l ∈ K0 , Az + Bzk Lkz −I + (∗) + 0 Pl + El Wl Ek H   "   ¯ ¯ ¯ ¯   0 G ¯ kz + B¯ zk L¯ kz −I + (∗) + −Pk + Ek Wk Ek < 0 k, l ∈ K1 , A 0 P¯l + E¯ l Wl E¯ k H¯

(63) (64)

(65) (66)

¯ H and H¯ being matrices of proper dimension. Let G = Gkl , G ¯ =G ¯ kl , H = Hkl and H¯ = H¯ kl for with G, G, k, l ∈ {1, 2, . . . , q} so the previous expressions render as:   Gkl Akz + Gkl Bzk Lkz + (∗) − Pk + EkT Wk Ek (∗) (67) < 0, k, l ∈ K0 , Hkl Akz + Hkl Bzk Lkz − GTkl −Hkl − HklT + Pl + ElT Wl El   ¯ kl B¯ zk L¯ kz + (∗) − P¯k + E¯ T Wk E¯ k ¯ kl A¯ kz + G G (∗) k (68) < 0, k, l ∈ K1 . ¯T −H¯ kl − H¯ klT + P¯l + E¯ lT Wl E¯ l H¯ kl A¯ kz + H¯ kl B¯ zk L¯ kz − G kl Therefore, LMIs (55) and (56) guarantee (67) and (68), respectively. Concerning “jumps” between regions of different kind, there are two cases: if k ∈ K1 and l ∈ K0 , then V < 0 if  T    −P¯k 0 x(t) ¯ x(t) ¯ V = < 0. (69) x(t ¯ + 1) x(t ¯ + 1) 0 Pˆl As in the previous case, we consider the following inequality from (29) and (28),   T  T   0 x(t) ¯ E¯ k Wk E¯ k x(t) ¯ ∀x(t), ¯ ∈ Ck , k ∈ K1 , ≥0 ∀x(t ¯ + 1), ∈ Cl , l ∈ K0 , x(t ¯ + 1) x(t ¯ + 1) 0 E¯ lT Wl E¯ l and the expression from (52)  "   k  x(t) ¯ =0 A¯ z + B¯ zk L¯ kz −I x(t ¯ + 1)

∀x(t), ¯ ∈ Ck , k ∈ K1 , ∀x(t ¯ + 1), ∈ Cl , l ∈ K1 .

(70)

(71)

¯ =G ¯ kl and H¯ = H¯ kl to obtain the By Finsler’s lemma, expressions (69), (70), and (71) are combined, choosing G inequality:   ¯ kl B¯ zk L¯ kz + (∗) − P¯k + E¯ T Wk E¯ k ¯ kl A¯ kz + G G (∗) k (72) < 0, k, l ∈ K1 . ¯T −H¯ kl − H¯ klT + Pˆl + E¯ lT Wl E¯ l H¯ kl A¯ kz + H¯ kl B¯ zk L¯ kz − G kl Hence, LMIs (57) guarantee (72). Similarly, if k ∈ K0 , l ∈ K1 , then V < 0 if T     x(t) ¯ x(t) ¯ −Pˆk 0 <0 V = x(t ¯ + 1) x(t ¯ + 1) 0 P¯l and considering the following inequality from (27) and (30)  T  T    x(t) ¯ x(t) ¯ 0 E¯ k Wk E¯ k ≥0 x(t ¯ + 1) x(t ¯ + 1) 0 E¯ lT Wl E¯ l

∀x(t), ¯ ∈ Ck , k ∈ K0 , ∀x(t ¯ + 1), ∈ Cl , l ∈ K1 ,

and the expression from (52) 1 The Finsler’s lemma [52] has been widely used in the TS context, especially in recent works. The interested reader is referred to [50,53].

(73)

(74)

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.19 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••



A¯ kz + B¯ zk Lˆ kz

−I





 " x(t) ¯ =0 x(t ¯ + 1)

∀x(t), ¯ ∈ Ck , k ∈ K0 , ∀x(t ¯ + 1), ∈ Cl , l ∈ K1 ,

we have  ¯ kl B¯ zk Lˆ kz + (∗) − Pˆk + E¯ T Wk E¯ k ¯ kl A¯ kz + G G k ¯T H¯ kl A¯ kz + H¯ kl B¯ zk Lˆ kz − G kl

19

(75)

 (∗) < 0, k, l ∈ K1 , −H¯ kl − H¯ klT + P¯l + E¯ lT Wl E¯ l

(76)

¯ =G ¯ kl , and H¯ = H¯ kl . Thus, LMIs (58) guarantee (76), concluding the proof. with G Example 7. Consider the following nonlinear model in [15]:     5 − x1 1 x1 x(t) + x(t + 1) = u(t) −1 −0.5 −2x1

(77)

Rewriting this model in a TS form within the compact set C = {x : |xi | ≤ β}, i ∈ {1, 2}, the following TS model is obtained: x(t + 1) = 

2 

hi (z(t)) (Ai x(t) + Bi u(t))

(78)

i=1

       β − x1 1 −β 1 β 5+β 5−β , A2 = , B1 = , B2 = , h1 (z(t)) = , h2 (z(t)) = 1 − −1 −0.5 −1 −0.5 2β −2β 2β h1 (z(t)). This model with the ordinary quadratic conditions, a solution is available when β < 1.36. Now consider the following 4-regions PWTS model obtained by dividing the compact subset C in the following subsets: with A1 =

C1 = {x : −β ≤ x1 ≤ −β/2, |x2 | ≤ β}, C3 = {x : 0 ≤ x1 ≤ β/2, |x2 | ≤ β},

C2 = {x : −β/2 ≤ x1 ≤ 0, |x2 | ≤ 2β},

C4 = {x : β/2 ≤ x1 ≤ β, |x2 | ≤ β}.

Therefore, the following 4-region PWTS model is obtained: x(t + 1) =

r 

  hki (z(t)) Aki x(t) + Bik u(t) ,

x(t) ∈ Ck ,

(79)

i=1

where zk = max(x1 ), zk = min (x1 ), hk1 = x∈Ck



Ak1

z k − x1

, hk = 1 − hk1 , and zk − zk 2     5 − zk 1 zk k k , B1 = , A2 = −1 −0.5 −2zk

x∈Ck

 1 zk = , −1 −0.5



B2k

 5 − zk = , −2zk

for k ∈ {1, . . . , 4}. The following matrices Ek , k ∈ {1, . . . , 4} were obtained using the PWLTOOL [39]         1 0 β −1 0 1 0 ¯ 1 0 −β/2 ¯ E4 = E1 = , E2 = , E3 = . −1 0 −β/2 −1 0 1 0 −1 0 β Conditions in Theorem 5 are expressed as BMIs and they have been solved using PENBMI toolbox in MATLAB [51]. This PWTS representation in (79) altogether with Theorem 5 allows establishing stability for β = 1.76, which outperforms the quadratic stability test for a single-compact TS representation. Thus, the following results were obtained:     L¯ 11 = −0.1240 0.2532 −0.0256 , L¯ 12 = −0.1284 0.2498 −0.0284 ,     L21 = −0.1143 0.0832 , L22 = −0.1142 0.0829 ,     L31 = −0.1520 −0.0760 , L32 = −0.1516 −0.0771 ,     L¯ 41 = −0.2558 −0.4118 0.0077 , L¯ 42 = −0.2559 −0.4120 0.0092 ,     2.6543 0.5597 2.4600 0.5223 , P3 = 1 × 103 , P2 = 1 × 103 0.5597 0.7896 0.5223 0.85681

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.20 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

20

Fig. 9. Results in Example 7.

⎡

⎤ 2.0624 0.0453 −0.2877 1.2241 −0.2101 ⎦ , P¯1 = 1 × 103 ⎣ 0.0453 −0.2877 −0.2101 −0.4524

⎡

⎤ 0.4200 0.1761 0.2029 P¯4 = 1 × 103 ⎣ 0.1761 1.7241 0.0048 ⎦ . 0.2029 0.0048 0.0217

Fig. 9(a) shows with solid lines the boundaries of each subcompact Ck , k ∈ {1, . . . , 4}; some Lyapunov curve levels q corresponding to (17) within C = k=1 Ck are also shown in dashed lines in the same figure. Four trajectories of the states from different initial conditions are shown to converge to the origin as it has been expected. In Fig. 9(b) it is shown the time evolution of the control law for the trajectories in Fig. 9(a). 4.3. Estimation of the domain of attraction In contrast with the stability analysis, the domain of attraction in local stabilization of nonlinear systems via TS models must take into account the following characteristics: 1. It can be affected by the control gains, which can be taken into account [27] or pre-fixed [37]. 2. It can be enlarged as long as the convex-sum property of the MFs hold, i.e., it is naturally bounded by the modeling area whereas for stability analysis this could be changed at convenience. Algorithm 1 as well as its modified version Algorithm 2 can be applied to the local stabilization problem, provided some slight modifications are made: Algorithm 3. Starting from the compact C0 with length α0 corresponding to the maximum s-hypercube on which the nonlinear system is quadratically stabilized with a single PDC control law, set i = 1 and α1 = α0 and perform the following steps: j

1. Consider s-hypercube compacts Ci of length αi indexed by j ∈ Ii and fully covering and surrounding Ci−1 by j sharing a side or a vertex with it, such that Ci ∩ Cik = ∅ for j = k. 2. Set Ci = Ci−1 . j 3. For each j ∈ Ii test whether Ci ∪ Ci under its corresponding PWTS representation renders Theorem 4 (if continj uous) or Theorem 5 (if discrete) feasible or not. If yes, redefine Ci = Ci ∪ Ci . Go to the next j until all the indexes in Ii have been run. Note that a piecewise control law of the form (44) is found at each iteration. j 4. If Ci includes all the regions Ci , j ∈ Ii (i.e., all of them altogether with Ci−1 are feasible in the PWTS sense), then set αi+1 as the length of Ci ; otherwise, set αi+1 as a half of the length of any Ci−1 . 5. The current estimation of the DA is given by the outermost Lyapunov set V (x) ≤ γ in Ci with V (x) as the PWLF (14) (continuous-time case) or (17) (discrete-time case). If this estimate is satisfactory or the LMIs are rendered

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.21 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

21

Fig. 10. Comparison: DA via Algorithm 3 (outermost black closed curve), DA via Th. 2 in [27] (blue closed curve). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

computationally intractable or the original modeling region is already reached, terminate the algorithm; otherwise, set i = i + 1 and go to step 1 for a finer estimation of the DA. Example 8. Consider the nonlinear model in [27,43]:     1 3 sin x1 + 1 −4 x(t) + u(t) x(t) ˙ = 4.5 sin x1 + 5.5 9.5 − 10.5 sin x1 −2

(80)

where the compact set C = R2 is the modeling area. Nonlinearity z1 = sin x1 is naturally bounded as z1 = −1 ≤ z1 ≤ z1 = 1 in R2 . Therefore, this model can be rewritten in the TS form (6) through the ordinary sector nonlinearity approach which corresponds to the methodology described above with q = 1 and         −2 −4 4 −4 1 1 , A2 = , B1 = , B2 = . A1 = 20 −2 −1 −2 1 10 Applying Algorithm 3 for as many iterations as to overcome results in [27] (12 iterations after the single PDC region is split in two around the origin), the estimation of the DA shown in Fig. 10 is obtained (outermost black closed curve); it clearly outperforms the estimation given in [27]. 5. Conclusion This paper introduced a novel stability analysis and controller design for nonlinear systems, both continuousand discrete-time, based on exact piecewise Takagi–Sugeno representations. It has been proved that the proposed methodology allows progressively better estimates of the domain of attraction to be found via easy-to-implement algorithms. An important open question to be considered in future research is: can the domain of attraction of an asymptotically stable equilibrium point x = 0 be asymptotically estimated with an arbitrary degree of accuracy via piecewise Lyapunov functions (and PWTS models) or is there an underlying limitation to this goal? Based on the results offered in this paper, it seems that these estimates are only subject to computational restrictions. A number of examples have been shown that easily overcome recent results on the subject. References [1] T. Takagi, M. Sugeno, Fuzzy identification of systems and its applications to modeling and control, IEEE Trans. Syst. Man Cybern. 15 (1) (1985) 116–132. [2] L. Wang, A Course in Fuzzy Systems and Control, Prentice–Hall, Inc., Upper Saddle River, NJ, USA, 1997. [3] T. Taniguchi, K. Tanaka, H. Wang, Model construction, rule reduction and robust compensation for generalized form of Takagi–Sugeno fuzzy systems, IEEE Trans. Fuzzy Syst. 9 (2) (2001) 525–537. [4] K. Tanaka, H. Wang, Fuzzy Control Systems Design and Analysis: A Linear Matrix Inequality Approach, John Wiley and Sons, New York, 2001.

JID:FSS AID:6943 /FLA

22

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.22 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

[5] Z. Lendek, T. Guerra, R. Babuska, B. De-Schutter, Stability Analysis and Nonlinear Observer Design Using Takagi–Sugeno Fuzzy Models, Springer-Verlag, Netherlands, 2010. [6] S. Boyd, L. Ghaoui, E. Feron, V. Belakrishnan, Linear Matrix Inequalities in System and Control Theory, Studies in Applied Mathematics, vol. 15, SIAM, Philadelphia, USA, 1994. [7] A. Sala, T. Guerra, R. Babuska, Perspectives of fuzzy systems and control, Fuzzy Sets Syst. 156 (3) (2005) 432–444. [8] G. Feng, A survey on analysis and design of model-based fuzzy control systems, IEEE Trans. Fuzzy Syst. 14 (5) (2006) 676–697. [9] H. Tuan, P. Apkarian, T. Narikiyo, Y. Yamamoto, Parameterized linear matrix inequality techniques in fuzzy control system design, IEEE Trans. Fuzzy Syst. 9 (2) (2001) 324–332. [10] X. Liu, Q. Zhang, New approaches to controller designs based on fuzzy observers for T–S fuzzy systems via LMI, Automatica 39 (9) (2003) 1571–1582. [11] A. Sala, C. Ariño, Asymptotically necessary and sufficient conditions for stability and performance in fuzzy control: applications of Polya’s theorem, Fuzzy Sets Syst. 158 (24) (2007) 2671–2686. [12] A. Kruszewski, A. Sala, T. Guerra, C. Arino, A triangulation approach to asymptotically exact conditions for fuzzy summations, IEEE Trans. Fuzzy Syst. 17 (5) (2009) 985–994. [13] Y. Blanco, W. Perruqueti, P. Borne, Stability and stabilization of nonlinear systems and Takagi–Sugeno fuzzy models, Math. Probl. Eng. 7 (3) (2001) 221–240. [14] K. Tanaka, T. Hori, H. Wang, A multiple Lyapunov function approach to stabilization of fuzzy control systems, IEEE Trans. Fuzzy Syst. 11 (4) (2003) 582–589. [15] T. Guerra, L. Vermeiren, LMI-based relaxed non-quadratic stabilization conditions for nonlinear systems in Takagi–Sugeno’s form, Automatica 40 (5) (2004) 823–829. [16] B. Rhee, S. Won, A new fuzzy Lyapunov function approach for a Takagi–Sugeno fuzzy control system design, Fuzzy Sets Syst. 157 (9) (2006) 1211–1228. [17] M. Bernal, T.M. Guerra, Generalized non-quadratic stability of continuous-time Takagi–Sugeno models, IEEE Trans. Fuzzy Syst. 18 (4) (2010) 815–822. [18] T. Guerra, M. Bernal, Strategies to exploit non-quadratic local stability analysis, Int. J. Fuzzy Syst. 14 (3) (2012) 372–379. [19] M. Johansson, A. Rantzer, K. Arzen, Piecewise quadratic stability of fuzzy systems, IEEE Trans. Fuzzy Syst. 7 (6) (1999) 713–722. [20] G. Feng, Stability analysis of discrete-time fuzzy dynamic systems based on piecewise Lyapunov functions, IEEE Trans. Fuzzy Syst. 12 (1) (2004) 22–28. [21] G. Feng, C. Chen, D. Sun, Y. Zhu, H∞ controller synthesis of fuzzy dynamic systems based on piecewise Lyapunov functions and bilinear matrix inequalities, IEEE Trans. Fuzzy Syst. 13 (1) (2005) 94–103. [22] M. Bernal, T.M. Guerra, A. Kruszewski, A membership-function-dependent approach for stability analysis and controller synthesis of Takagi– Sugeno models, Fuzzy Sets Syst. 160 (19) (2009) 2776–2795. [23] K. Tanaka, H. Ohtake, H. Wang, A descriptor system approach to fuzzy control system design via fuzzy Lyapunov functions, IEEE Trans. Fuzzy Syst. 15 (3) (2007) 333–341. [24] K. Guelton, T. Bouarar, N. Manamanni, Robust dynamic output feedback fuzzy Lyapunov stabilization of Takagi–Sugeno systems – a descriptor redundancy approach, Fuzzy Sets Syst. 160 (19) (2009) 2796–2811. [25] M. Bernal, A. Sala, A. Jaadari, T.M. Guerra, Stability analysis of polynomial fuzzy models via polynomial fuzzy Lyapunov functions, Fuzzy Sets Syst. 185 (1) (2011) 5–14. [26] J. Pitarch, A. Sala, C. Ariño, F. Bedate, Estimación del dominio de atracción de sistemas no lineales mediante modelos borrosos polinomiales (in Spanish), Rev. Iberoam. Autom. Inf. Ind. 9 (2) (2012) 152–161. [27] D. Lee, D. Kim, Relaxed LMI conditions for local stability and local stabilization of continuous-time Takagi–Sugeno fuzzy systems, IEEE Trans. Cybern. 44 (3) (2014) 394–405. [28] T. Hori, K. Tanaka, H. Wang, A piecewise Takagi–Sugeno fuzzy model construction and relaxation of stability conditions, in: Proceedings of the 41st IEEE Conference on Decision and Control, vol. 2, Las Vegas, USA, 2002, pp. 2149–2150. [29] H. Ohtake, K. Tanaka, H. Wang, Piecewise nonlinear control, in: Proceedings of the 42nd IEEE Conference on Decision and Control, vol. 5, Maiu, USA, 2003, pp. 4735–4740. [30] V. Campos, F. Souza, L. Torres, R. Palhares, New stability conditions based on piecewise fuzzy Lyapunov functions and tensor product transformations, IEEE Trans. Fuzzy Syst. 21 (4) (2013) 748–760. [31] T. Gonzalez, M. Bernal, R. Marquez, Stability analysis of nonlinear models via exact piecewise Takagi–Sugeno models, in: Proceedings of the 19th IFAC World Congress, Cape Town, South Africa, 2014, pp. 1–6. [32] T. Guerra, A. Kruszewski, J. Lauber, Discrete Tagaki–Sugeno models for control: where are we?, Annu. Rev. Control 33 (1) (2009) 37–47. [33] T. Guerra, A. Kruszewski, M. Bernal, Control law proposition for the stabilization of discrete Takagi–Sugeno models, IEEE Trans. Fuzzy Syst. 17 (3) (2009) 724–731. [34] T. Guerra, M. Bernal, K. Guelton, S. Labiod, Nonquadratic local stabilization for continuous-time Takagi–Sugeno models, Fuzzy Sets Syst. 201 (2012) 40–54. [35] J. Pitarch, A. Sala, C. Ariño, Closed-form estimates of the domain of attraction for nonlinear systems via fuzzy-polynomial models, IEEE Trans. Cybern. 44 (4) (2014) 526–538. [36] H. Zhang, G. Feng, Stability analysis and H∞ controller design of discrete-time fuzzy large-scale systems based on piecewise Lyapunov functions, IEEE Trans. Syst. Man Cybern., Part B, Cybern. 38 (5) (2008) 1390–1401. [37] C. Arino, E. Perez, A. Sala, F. Bedate, Polytopic invariant and contractive sets for closed-loop discrete fuzzy systems, J. Franklin Inst. 351 (7) (2014) 3559–3576. [38] J. Qiu, H. Tian, Q. Lu, H. Gao, Nonsynchronized robust filtering design for continuous-time T–S fuzzy affine dynamic systems based on piecewise Lyapunov functions, IEEE Trans. Cybern. 43 (6) (2013) 1755–1766.

JID:FSS AID:6943 /FLA

[m3SC+; v1.218; Prn:19/11/2015; 10:04] P.23 (1-23)

T. González, M. Bernal / Fuzzy Sets and Systems ••• (••••) •••–•••

23

[39] S. Hedlund, M. Johansson, Pwltlool: a matlab toolbox for analysis of piecewise linear systems, 1999. [40] G. Feng, Controller synthesis of fuzzy dynamic systems based on piecewise Lyapunov functions, IEEE Trans. Fuzzy Syst. 11 (5) (2003) 605–612. [41] F. Shirani, M. Yazdanpanah, B. Araabi, Comments on “controller synthesis of fuzzy-dynamic systems based on piecewise Lyapunov function”, IEEE Trans. Fuzzy Syst. 18 (1) (2010) 227–228. [42] H. Khalil, Nonlinear Systems, Pearson Education, New Jersey, USA, 2000. [43] D. Lee, J. Park, Y. Joo, A fuzzy Lyapunov function approach to estimating the domain of attraction for continuous-time Takagi–Sugeno fuzzy systems, Inf. Sci. 185 (1) (2012) 230–248. [44] A. Chakraborty, P. Seiler, G. Balas, Nonlinear region of attraction analysis for flight control verification and validation, Control Eng. Pract. 19 (4) (2011) 335–345. [45] G. Chesi, Domain of Attraction Analysis and Control via SOS Programming, Springer, Berlin, Germany, 2011. [46] E. Tognetti, R. Oliveira, P. Peres, Selective H2 and H∞ stabilization of Takagi–Sugeno fuzzy systems, IEEE Trans. Fuzzy Syst. 19 (5) (2011) 890–900. [47] H. Zhang, X. Xie, Relaxed stability conditions for continuous-time T–S fuzzy-control systems via augmented multi-indexed matrix approach, IEEE Trans. Fuzzy Syst. 19 (3) (2011) 478–492. [48] G. Feng, G. Lu, Design of discrete time fuzzy controller based on piecewise Lyapunov functions, in: Proceedings of the 4th World Congress on Intelligent Control and Automation, vol. 4, Shanghai, China, 2002, pp. 3247–3251. [49] C. Chen, G. Feng, D. Sun, X. Guan, H∞ output feedback control of discrete-time fuzzy systems with application to chaos control, IEEE Trans. Fuzzy Syst. 13 (4) (2005) 531–543. [50] J. Pan, T. Guerra, S. Fei, A. Jaadari, Nonquadratic stabilization of continuous T–S fuzzy models: LMI solution for a local approach, IEEE Trans. Fuzzy Syst. 20 (3) (2012) 594–602. [51] M. Koˇcvara, M. Stingl, Penbmi user’s guide. version 2.1, Software manual, PENOPT GbR, 2006. [52] P. Finsler, Über das Vorkommen definiter und semi-definiter Formen in scharen quadratischer Formen, Comment. Math. Helv. 9 (1) (1937) 188–192. [53] A. Jaadari, T. Guerra, A. Sala, M. Bernal, Finsler’s relaxation for local H-infinity controller design of continuous-time Takagi–Sugeno models via non-quadratic Lyapunov functions, in: Proceedings of the 2013 American Control Conference, Washington, DC, USA, 2013, pp. 5668–5673.