Constructive proof of Lagrange stability and sufficient – Necessary conditions of Lyapunov stability for Yang–Chen chaotic system

Applied Mathematics and Computation 309 (2017) 205–221

Contents lists available at ScienceDirect

Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc

Constructive proof of Lagrange stability and sufficient – Necessary conditions of Lyapunov stability for Yang–Chen chaotic systemR Xiaoxin Liao a,b, Guopeng Zhou a,∗, Qigui Yang c, Yuli Fu d, Guanrong Chen e a

Institute of Engineering and Technology, Hubei University of Science and Technology, Xianning 437100, PR China College of Automation, Huazhong University of Science and Technology, Wuhan 430074, PR China c School of Mathematics, South China University of Technology, Guangzhou 510641, PR China d School of Electronics and Information Engineering, South China University of Technology, Guangzhou 510641, PR China e Department of Electronic Engineering, City University of HongKong, Hong Kong Special Administrative Region, PR China b

a r t i c l e

i n f o

Keywords: Yang–Chen system Lagrange stability Global exponential attractive set Lyapunov stability Branch

a b s t r a c t This paper studies the stability problem of Yang–Chen system. By introducing different radial unbounded Lyapunov functions in different regions, global exponential attractive set of Yang–Chen chaotic system is constructed with geometrical and algebraic methods. Then, simple algebraic sufficient and necessary conditions of global exponential stability, global asymptotic stability, and exponential instability of equilibrium are proposed. And the relevant expression of corresponding parameters for local exponential stability, local asymptotic stability, exponential instability of equilibria are obtained. Furthermore, the branch problem of the system is discussed, some branch expressions are given for the parameters of the system. © 2017 Elsevier Inc. All rights reserved.

1. Introduction Since Lorenz chaotic system was found in 1963 [1,2], chaotic systems have attracted many researchers’ interests, and great progress has been made in many areas. Not only basic properties and dynamical behaviors of Lorenz system have been widely studied, but also chaotic control, chaotic synchronization, and other applications have been taken cared of [3– 28]. Furthermore, some new chaotic systems which are not equivalent to Lorenz system were found, i.e., Rössler system [4], Chua system [5], Chen system [8,9], Lü–Chen system [10], Yang–Chen system [13], etc. The authors in [31–36] study the chaotic systems by using FPGA method to realize the simulation. It is generally believed that the conditions for producing chaos of a continuous dynamical system include: (I) the system has at least one positive Lyapunov exponent in a small area; (II) the system is ultimately bounded (or Lagrange stable), namely, the trajectory of the system far away from the equilibrium point converges to a specific bounded set. Furthermore, computing Lyapunov exponent is significant on the condition that (II) holds. Therefore, the proof of ultimate boundedness for dynamical systems is a core and critical problem. Soviet academician Leonov is the first one who studied this problem. R This work is partially supported by the Joint Key Grant of National Natural Science Foundation of China and Zhejiang Province (U1509217), Hubei Province Science and Technology Support Program (2015BAA001), Hubei SMEs Innovation Fund Project (2015DAL069). ∗ Corresponding author. E-mail addresses: [email protected], [email protected] (G. Zhou).

http://dx.doi.org/10.1016/j.amc.2017.03.033 0 096-30 03/© 2017 Elsevier Inc. All rights reserved.

206

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

Fig. 1. The phase diagram of Yang–Chen system.

He respectively obtained a cylindrical estimator and one oval estimator for Lorenz system [14–16,25]. Based on the above results, we not only gave a simplified proof for of ultimate boundedness of Lorenz system, but also improved and extended the well-known results [23,24]. Particularly, we first proposed a new concept of global exponential stable in sense of Lagrange [23]. Recently, by using geometric and algebraic methods, we proposed a constructive proof of globally exponentially attractive set for Chen system [18]. Lyapunov stability of equilbria for nonlinear dynamical systems is important. In the literature [25], we have completed Lyapunov stability analysis of equilibrium point of Lorenz system. Some very simple algebraic sufficient and necessary conditions of global exponential stability, global asymptotical stability, and instability were given, and the results were applied to chaotic control very well. But for Yang–Chen system, does it has the same result? In this paper, we focus on studying Lagrange stability and Lyapunov stability of Yang–Chen system. For convenience of our discussion, Yang–Chen system is described as follows [13]:

⎧ dx ⎪ = a ( y − x ), ⎪ ⎪ dt ⎪ ⎪ ⎨ dy

⎪ dt ⎪ ⎪ ⎪ ⎪ d ⎩ z dt

= dx − xz,

(1)

= xy − bz,

where x, y, z ∈ R are state variables, a, b, d are parameter constants, and a, b are always assumed to be positive values. Let X (t )√= (x(t√ ), y(t ), z(t )) with X (t ) = X (t, t0 , X0 ), and X0 is the initial value of X(t). It has three equilibria S0 (0, 0, 0 ), S± (± bd, ± bd, d ). When a = 10, b = 1, d = 16, system (1) has a chaotic attractor as it is shown in Fig. 1. Since it has the property of invariance and symmetry (i.e., (−x, −y, z ) → (x, y, z )), asymptotic behavior of Yang–Chen system will be studied by transforming three-dimensional state (x, y, z) into two-dimensional plane (y, z). By introducing different radial unbounded Lyapunov functions in different regions, global exponential attractive set of Yang–Chen chaotic system is constructed with geometrical and algebraic methods. Then, simple algebraic sufficient and necessary conditions of global exponential stability, global asymptotic stability, and instability of equilibrium S0 are proposed. And we also obtain the relevant expression of the corresponding parameters for local exponential stability, local asymptotic stability, exponential instability of equilibria S± . Furthermore, the branch problem is discussed, some branch expressions are given for the parameters of Yang–Chen system. Lagrange stability and Lyapunov stability conditions presented here are very useful for studying how chaos produced, chaos control, chaos synchronization, and other applications. The rest of this paper is organized as follows. In the second section, some necessary preparations are given. In Section 3, global exponential attractive set will be presented with a constructive proof. Algebraic sufficient and necessary conditions of Lyapunov stability of equilibrium S0 are proposed in Section 4. In Section 5, the conditions of local stability of equilibria S± are presented. In Section 6, the branch value expressions of different parameters are obtained. Section 7 is the simulation to show the correctness of the results presented in the previous sections. And the final section is the conclusion. 2. Preliminaries We will discuss the properties of Lagrange stability, Lyapunov stability and bifurcation for system (1). For this purpose, some definitions and lemmas are presented in this section. Definition 1. If there exists a compact set  ⊂ R3 for all X0 ∈ c ⊆R3 and ρ (X (t ), )  limt→+∞ infY ∈ X (t ) − Y = 0 holds, then the compact set  is called a global attractive set of system (1). That is, all solutions of system (1) are ultimately bounded. The complementary set c of  is called global attractive domain.

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

207

System (1) is also called global asymptotic stability in sense of Lagrange or ultimate boundedness in sense of dissipation if it satisfies the conditions of Definition 1. Furthermore, if ∀X0 ∈ 0 ⊆ ⊂ R3 , ∀t ≥ t0 , X(t, t0 , X0 ) ∈ 0 , then, 0 is called positive invariant set. Definition 2. If there exist a compact set  ⊂ R3 and two constants M(X0 ) > 0, α > 0, ∀X0 ∈ c , it holds that ρ (X (t ), ) ≤ M (X0 )e−α (t−t0 ) . Then, system (1) is globally exponentially attractive and  is called globally exponentially attractive set. Since the existence of  in Definition 2 is difficult to be verified, and Lyapunov function method is a powerful tool to study asymptotic property of nonlinear dynamic systems, another definition which is easier to be proved than Definition 2 is introduced in the following. Definition 3. If there is a generalized positive definite and radially unbounded Lyapunov function V(X(t)) and two positive constants L > 0, α > 0. ∀t ≥ t0 , it holds that 0 < V (X (t )) − L ≤ (V (X0 ) − L )e−α (t−t0 ) . Then, system (1) is said to be globally exponentially attractive or globally exponentially stable in sense of Lagrange,  = {X |V (X (t )) ≤ L, t ≥ t0 } is globally exponentially attractive set. Lemma 1. If the trajectory y(t) of system (1) has a globally exponentially attractive set, that is, there exist two constants k > 0 and α > 0, such that the following estimation holds:

|y(t )| ≤ k + |y0 |e−α (t−t0 ) , ∀ t ≥ t0 .

(2)

Then, the trajectory x(t) of system (1) also has a global exponential attractive set, and it has the following estimation:

|x(t )| ≤

αk α¯ + |x0 |e−α (t−t0 ) + a|y0 |e− 2 (t−t0 ) , ∀t ≥ T0 , α¯

(3)

where α¯ = min{a, α}, T0 ≥ t0 . Proof. From the first equation of system (1), one obtains that

x(t ) = x0 e−a(t−t0 ) +



t

t0

e−a(t−τ ) ay(τ )dτ .

Consider inequality (2), one has

|x(t )| ≤ |x0 |e−a(t−t0 ) +



t

t0

ae−α¯ (t−τ ) (k + |y0 |e−α¯ (τ −t0 ) )dτ

ak −α¯ t α¯ t = |x0 |e−a(t−t0 ) + e (e − eα¯ t0 ) α¯  t +a|y0 | e−α¯ (t−t0 ) dτ t0

= |x0 |e−α (t−t0 ) + α¯

ak (1 − e−α¯ (t−t0 ) ) α¯ α¯

+a|y0 |e− 2 (t−t0 ) e− 2 (t−t0 ) (t − t0 ). α¯

α¯

Since limt→+∞ e− 2 (t−t0 ) (t − t0 ) = 0, there exists a constant T0 ≥ t0 , ∀t ≥ T0 , it holds that e− 2 (t−t0 ) (t − t0 ) ≤ 1. Then, it has

|x(t )| ≤

ak α¯ + |x0 |e−α (t−t0 ) + a|y0 |e− 2 (t−t0 ) , α¯

The proof is complete.

∀t ≥ T0 .



Lemma 2. If x(t) and y(t) have globally exponentially attractive sets, respectively. That is, they have the following estimations:

|x(t )| ≤ k + |x0 |e−α (t−t0 ) , ∀t ≥ t0 , |y(t )| ≤ k + |y0 |e−α (t−t0 ) , ∀t ≥ t0 , where k > 0, α > 0 are constants. Then, the solution z(t) of system (1) also has a globally exponentially attractive set, that is, it has the following estimation.

|z(t )| ≤

k2 + |z0 |e−b(t−t0 ) b¯ b¯

+(k|x0 | + k|y0 | + |x0 ||y0 | )e− 2 (t−t0 ) , ∀t ≥ T1 where b¯ = min{b, α}.

(4)

208

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

Proof. From the third equation of system (1), one has

z(t ) = z0 e−b(t−t0 ) +



t

t0

Hence,

|z(t )| ≤ |z0 |e−b(t−t0 ) +

e−b(t−τ ) x(τ )y(τ )dτ . 

t

t0

e−b(t−τ ) [k + |x0 |e−α (τ −t0 ) ]

[k + |y(x0 )|e−α (τ −t0 ) ]dτ  t = |z0 |e−b(t−t0 ) + e−b(t−τ ) t0

[k + k(|x0 | + |y0 | )e−α (τ −t0 ) + |x0 y0 |e−2α (τ −t0 ) ]dτ 2

k2 (1 − e−b(t−t0 ) ) b  t ¯ +[k|x0 | + k|y0 | + |x0 y0 |] e−b(t−t0 ) dτ

≤ |z0 |e−b(t−t0 ) +

t0

= |z0 |e−b(t−t0 ) +

k2 (1 − e−b(t−t0 ) ) b

+[k|x0 | + k|y0 | + |x0 y0 |]e−b(t−t0 ) (t − t0 ). ¯

b¯

Similarly, there exists a constant T1 ≥ t0 , ∀t ≥ T1 , such that e− 2 (t−t0 ) (t − t0 ) ≤ 1. Then, it holds that

|z(t )| ≤ |z0 |e−b(t−t0 ) +

k2 b b¯

+[k|x0 | + k|y0 | + |x0 ||y0 |]e− 2 (t−t0 ) , t ≥ T1 . 

The proof is complete.

Remark 1. From these two lemmas, one can see that if y(t) is globally exponentially stable in sense of Lagrange, then, x(t) is globally exponentially stable in sense of Lagrange. If x(t) and y(t) are globally exponentially stable in sense of Lagrange, then, z(t) is also globally exponentially stable in sense of Lagrange. In the process of the following main theorem, let |x| > c = const, it can be proved that (y(t), z(t)) is globally exponentially attractive. However, global exponential attractiveness of y(t) implies global exponential attraction of |x(t)|. Then, (x(t), y(t), z(t)) is globally exponentially stable in sense of Lagrange. 3. A constructive proof of globally exponentially attractive sets In this section, we will focus on globally exponentially attractive sets of system (1), a constructive proof will be proposed in the following process. At the outset, some important matrices are given as follows.

⎡

P1 

⎢ 1 ⎣− 2

P2 

P3 

1

0

⎤

1 2

⎡ ξ − 2a − d ⎥ 0 0 ⎦ > 0, Q1  ⎣

0

1 2

− 12

1+ε 0

0

1 − ε¯

0

0

1 2

1 2

0



> 0, Q2 



0

−1 0

0 − 2b



0

0

−a

0

0

− 2b

⎤

⎦ < 0,

< 0,

> 0, Q3  Q2 ,

where ε , ε¯ ∈ (0, 1 ) and ξ ∈ (0, 2a + d ) are constants selected appropriately such that P2 , P3 are symmetric positive definite, Q1 is negative definite. λM (Pi ), λM (Qi ), i = 1, 2, 3, are maximum eigenvalues of the corresponding matrices, respectively; λm (Pi ), i = 1, 2, 3, is the minimum eigenvalue of Pi √ . √ Let (y, z ) = (0, 3a + d ) be the center, and R1 =

λM (P1 ) 2b(3a+d ) ¯ √ , R1 = R1 + ε be the radii, two circles on (y, z) plane can 2 | λM ( Q 1 ) |

be described as follows:

⎧ 1 ⎪ ⎨ C : [y2 + (z − (3a + d ))2 ] = R21 , 2

⎪ ⎩C ε : 1 [y2 + (z − (3a + d ))2 ] = R¯ 21 . 2

(5)

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

209

These circles are respectively in quadrant I, II, III, IV of (y, z) plane, which intersect with the line z = η and z = ξ at

A1 (y1 , z1 ), A2 (y2 , z2 ), Aε1 (yε1 , z1ε ), Aε2 (yε2 , z2ε ). The equations of two lines connecting OA1 and OAε1 are z = ly, z = l¯y (l¯ < l ), respectively. And the equations of two lines connecting OA2 and OAε are z = −hy, z = −h¯ y (h¯ < h, 0 < η ξ ), respectively. 2

Draw two ellipses through two points A1 (y1 , z1 ) and Aε1 (yε1 , z1ε ), respectively. The corresponding equations are given as follows:

⎧ 1+ε 2 1 2 1+ε 2 1 2 ⎪ y + z = y1 + z1  R22 , ⎨O1 : 2

2

2

2

(6)

⎪ ⎩Oε : 1 + ε y2 + 1 z2 = 1 + ε (yε )2 + 1 (zε )2  R¯ 2 . 2 1 1 1 2

2

2

2

where R¯ 2 = R2 + ε . O1 and Oε1 intersect with axis z at A3 (0, z3 ), Aε3 (0, z3ε ), respectively. Similarly, draw two ellipses through two points A2 (y2 , z2 ) and Aε2 (yε2 , z2ε ), respectively. And the corresponding equations are given by:

⎧ 1 − ε¯ 2 1 2 1 − ε¯ 2 1 2 ⎪ y + z = y2 + z2  R23 , ⎨O2 : 2

2

2

2

(7)

⎪ ⎩Oε : 1 − ε¯ y2 + 1 z2 = 1 − ε¯ (yε )2 + 1 (zε )2  R¯ 2 . 3 2 2 2 2

2

2

2

where R¯ 3 = R3 + ε . O2 and Oε2 intersect with axis z at A4 (0, z4 ), Aε4 (0, z4ε ), respectively. Choose

⎧ 1 ⎪ 1 = {(y, z )| y2 + ⎪ ⎪ 2 ⎪ ⎪ ⎪ ⎪ 1 ⎪ ⎪ ⎨2 = {(y, z )| y2 + 2

1 (z − 3a − d )2 ≤ R21 , y ≥ 0, z ≤ η}, 2 1 (z − 3a − d )2 ≤ R21 , y < 0, z ≤ ξ }, 2

⎪ 1+ε 2 ⎪ ⎪ 3 = {(y, z )| y + ⎪ ⎪ 2 ⎪ ⎪ ⎪ ⎪ ⎩ = {(y, z )| 1 − ε¯ y2 + 4 2

and

⎧ 1 ⎪ ε1 = {(y, z )| y2 + ⎪ ⎪ 2 ⎪ ⎪ ⎪ ⎪ 1 ⎪ ⎪ ⎨ε2 = {(y, z )| y2 + 2

1 2 z ≤ R22 , y ≥ 0, z ≥ ly}, 2 1 2 z ≤ R23 , y < 0, z ≥ −hy}, 2

1 (z − 3a − d )2 ≤ R¯ 21 , y ≥ 0, z ≤ η}, 2 1 (z − 3a − d )2 ≤ R¯ 21 , y < 0, z ≤ ξ }, 2

⎪ 1+ε 2 ⎪ ⎪ ε3 = {(y, z )| y + ⎪ ⎪ 2 ⎪ ⎪ ⎪ ⎪ ⎩ε = {(y, z )| 1 − ε¯ y2 + 4 2

1 2 z ≤ R¯ 22 , y ≥ 0, z ≥ l¯y}, 2 1 2 z ≤ R¯ 23 , y < 0, z ≥ −h¯ y}. 2

Let  be the closed area of

1







and ε be the closed area of 



(see in Fig. 2). Then,  =



3 4 and ε = ε1 ε2 ε3 ε4 . ε 1−ε¯ 1+ε 1−ε¯ c c 2 Let G = {x||x| ≤ max{ 1+ εl d , εh d }}, G = {x||x| > max{ εl d , εh d }}, ε = R − ε . Then, we have the following result. 2

Theorem 1. Gc × cε is a global attractive domain of system (1), G × ε is a globally exponentially attractive set and positive invariant set of system (1). Thus, system (1) is an ultimately bounded chaotic system. Proof. The proof of Theorem 1 is divided into four steps. Step 1: When z ≤ ξ , construct a positive definite and radially unbounded Lyapunov function as follows:

V1 = x2 +

 =

1 2 1 y − xy + (z − 3a − d )2 2 2

x y z − 3a − d

T  P1



x y . z − 3a − d

(8)

210

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

Fig. 2. The schematic diagram of compact set .

Take the derivative of V1 along system (1), one has

dV1 dx dy dx dy dz = 2x +y −y −x + ( z − 3a − d ) dt dt dt dt dt dt −(2a + d )x2 − ay2 + x2 z − bz2 + (3a + d )bz ≤ −(2a + d − ξ )x2 − ay2 + x2 (z − ξ ) b b − ( z − 3a − d )2 + ( 3a + d )2 2 2



≤

x y z − 3a − d

T



Q1





x b y + ( 3a + d )2 2 z − 3a − d

x y ≤ λM (Q1 ) z − 3a − d

T 



x b y + ( 3a + d )2 2 z − 3a − d

λM (Q1 ) λM (P1 ) b [V − ( 3a + d )2 ] λM (P1 ) 1 |λM (Q1 )| 2 λM (Q1 ) ≤ [V − R21 ]. λM (P1 ) 1 ≤

Hence,

 ⎧ dV1  ⎪ ≤ ⎪ ⎨ dt  V1 =R21 dV1  ⎪ ⎪ ≤ ⎩ dt  2 V1 =R¯ 1

λM (Q1 ) 2 [R − R21 ] = 0, λM (P1 ) 1 λM (Q1 ) ¯ 2 [R − R21 ] < 0. λM (P1 ) 1

Thus, when z ≤ ξ , the trajectory X(t) is exponentially decreasing and entering into V1 ≤ R¯ 21 from V1 > R¯ 21 , and it will stay in the area of V1 < R¯ 21 . Step 2: Because of the symmetry property of (x, y) of system (1), namely, let (x, y ) = (−x, −y ), system (1) is invariant, we just consider the trend on x ≥ 0. Consider the positive definite and radially unbounded Lyapunov function on Quadrant I of (y, z) plane with z ≥ η:

V2 = =

1+ε 2 1 y + ( z − ( 1 + ε )d )2 2 2

y z − ( 1 + ε )d

T P2



y . z − ( 1 + ε )d

(9)

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

Consider x as a parameter with x ≥ (1), one has

1

ε l¯

211

¯

− 2bεl = x∗ , take the derivative of V2 along the second and third equations of system

dV2 =(1 + ε )yy˙ + (z − (1 + ε )d )z˙ dt b b b = − ε xyz − (z − (1 + ε )d )2 + ((1 + ε )d )2 − z2 2 2 2 b b ≤ − ε xl¯y2 − (z − (1 + ε )d )2 + ((1 + ε )d )2 2 2 ⎡ ⎤ b

T

 −ε xl¯ − l¯2 0 y y 2 ⎣ ⎦ = b z − ( 1 + ε )d z − ( 1 + ε )d 0 2 b + ((1 + ε )d )2 2

y ≤ z − ( 1 + ε )d

T

Q2

≤λM (Q2 )



b y + ((1 + ε )d )2 z − ( 1 + ε )d 2

y z − ( 1 + ε )d

T

y z − ( 1 + ε )d



b ((1 + ε )d )2 2 λM (Q2 ) λM (P2 ) b ≤ [V − ((1 + ε )d )2 ]. λM (P2 ) 2 |λM (Q2 )| 2 +

Hence,

V2 (y, z ) −

λM (P2 ) b ((1 + ε )d )2 ≤ |λM (Q2 )| 2

 λM (Q2 ) (t − t0 ) λM (P2 ) b 2 V2 (y0 , z0 ) − ((1 + ε )d ) e λM (P2 ) . |λM (Q2 )| 2



Then, the trajectory (y, z) is exponentially decreasing and enters into  from the outside through

.

Step 3: We construct a generalized positive definite and radially unbounded Lyapunov function for y < 0, z > ξ > 0 in Quadrant II of (y, z) plane.

V3 =

1 − ε¯ 2 1 y + (z − (1 − ε¯ )d )2 2 2

=

y z − (1 − ε¯ )d

T P3



y , z − (1 − ε¯ )d

where ε¯ ∈ (0, 1 ). The derivative of V3 along system (1) is given by:

b b dV3 b =ε¯ xyz − (z − (1 − ε¯ )d )2 + ((1 − ε¯ )d )2 − z2 dt 2 2 2 b b b ≤ − ε¯ xh¯ y2 − (z − (1 − ε¯ )d )2 + ((1 − ε¯ )d )2 − z2 2 2 2 ⎡ ⎤

T

 ¯ − b h¯ 2 h − ε ¯ x 0 y y 2 ⎣ ⎦ ≤ b z − (1 − ε¯ )d z − (1 − ε¯ )d 0 − 2 b + ((1 − ε¯ )d )2 2

≤

y z − (1 − ε¯ )d

T

Q3



b y + ((1 − ε¯ )d )2 . z − (1 − ε¯ )d 2

(10)

212

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

Hence, when x ≥

V3 (y, z ) −

1

ε h¯

− 2bε¯ h¯  x∗ , it has

λM (P3 ) b ((1 − ε¯ )d )2 |λM (Q3 )| 2

 − λM (P2 ) (t − t ) 0 λM (P3 ) b 2 ≤ V3 (y0 , z0 ) − ((1 − ε¯ )d ) e |λM (Q2 )| . |λM (Q3 )| 2 

Thus, (y, z) will exponentially decrease and enter into  from the outside of

.

Step 4: Adjust the following parameters appropriately, one can make z4 ≥ z3 , which makes point A4 (0, z4 ) above A3 (0, z3 ) in (y, z) plane. (a) Increase ξ close to 2a + d, h can be increased, so can z4 . (b) Decrease η close to zero such that l can be decreased to zero, thus z3 can be decreased. (c) Making appropriate adjustment for ε and ε¯ can make z4 ≥ z3 . Therefore, we might as well assume z4 ≥ z3 . If z4 = z3 , the result is obtained. If z4 > z3 , we only need prove that the trajectory (y(t), z(t)) of system (1), which is close to line A3 A4 with y(t ) = 0, will pass through the line A3 A4 from right to left. Hence, (y(t), z(t)) always trends from external of  to inside. In fact, plug y = 0 into the second equation of system (1), it has

dy = dx − xz = x(d − z ) < 0. dt

(11)

Since z > 0, the trajectory (y(t), z(t)) on line A3 A4 trends from right to left. By the Lefeshtz principle of continuity and extension [29], the trajectory (y(t), z(t)) on line A3 A4 also decreases exponentially. The proof is complete.  4. Simple algebraic sufficient and necessary condition and application for Lyapunov stability of equilibrium S0 There are three equilibria of system (1), it is given by:





S0  (0, 0, 0 ), S±  (± bd, ± bd, d ). The corresponding linearized coefficient matrices of system (1) at equilibria S0 , S± are respectively described as follows:



A0 

−a d 0



a 0 0

0 0 , A±  −b



−a 0 √ ± bd

a 0 √ ± bd



0 √ ∓ bd . −b

Theorem 2. Sufficient and necessary condition for global exponential stability of S0 in sense of Lyapunov is d < 0. Proof. Sufficiency. The proof process contains three parts. Part 1. Since a > 0, d < 0, one can choose some positive numbers α > β > 0 such that the following conditions hold. Condition I: α > β 2 ; Condition II: α + β = − da ; Condition III: aα + 2dβ > 0. These positive numbers α , β are easily obtained. If let

α=−

10n − 1 d 1 d ,β = − n , 10n a 10 a

one has

  10n − 1 1 d d α+β = − − =− . n n 10

10

a

a

By choosing n appropriately, one has

β 2 (− 101n da )2 −1 d = < 1, n (10n − 1 ) a d α − 1010n −1 10 n a 10n − 1 d 1 d − 2d n 10n a 10 a (10n − 1 )a + 2d d − > 0. 10n a aα + 2dβ − a

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

213

Fig. 3. The schematic diagram of rectangular area ABDEA.

Consider the following Lyapunov function candidate:

V4 =

α 2

x2 +

 T x y z

=

y2 z2 − β xy + 2 2

⎡

α

−

⎢ 2 ⎢ β ⎢− ⎣ 2

β

2 1 2 0

0

⎤    T   ⎥ x x x ⎥ 0 ⎥ y  y P4 y . ⎦ z z z 0

(12)

1 2

When z < −d, take the derivative of V4 along system (1), it yields

dV4 = α xx˙ + yy˙ − β x˙ y − β xy˙ + zz˙ dt ≤ −(aα + 2dβ )x2 − aβ y2 − bz2 + β x2 (z + d )

⎡ ⎤T ⎡ x

≤ ⎣y⎦ ⎣

−(aα + 2dβ )

0

0

−aβ

0

0

z

⎡ ⎤T

⎡ ⎤

x



0

⎤⎡ ⎤ x

0 ⎦⎣y⎦ −b

z

x

⎣y⎦ Q4 ⎣y⎦ < 0. z

z

Then, one has λM (Q4 )

V4 (X ) ≤ V4 (X0 )e λM (P4 )

(t−t0 )

.

(13)

Therefore, when z < −d, (x, y) will exponentially converge to (0, 0) or into the area of z ≥ −d. In Theorem 1, we proved that system (1) has a global exponential attraction . Now, construct a square ABCDEF in (y, z) plane, and establish a coordinate transformation yz  which is a one-one mapping from  to square ABCDEF (see in Fig. 3). Draw line z = η which intersects line AB at point A1 (y1 , η); and draw line z = |d| which intersects line DE at point A¯ 1 (y¯ 1 , |d| ). The linear equations of OA1 and OA¯ 1 are z = ly and z = −hy, respectively. Part 2. On the area of x ≥ 0, y ≥ 0, z ≥ ly (in OA1 BCO), construct a positive definite Lyapunov function:

V5 =

1+ε 1+ε 2 1 2 |d |x2 + y + z 2a 2 2

 T  1+ε =

x y z

2a

|d |

0 0

0

1+ε 2

0

0 0 1 2

  x y z

 T   

x y z

x P5 y . z

(14)

214

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

When x ≥ x∗ > 0 with x∗ a positive constant, take the derivative of V5 along system (1), it holds that

dV5 = (1 + ε )|d|xx˙ + (1 + ε )yy˙ + zz˙ dt ≤ −(1 + ε )|d|x2 + (1 + ε )|d|xy + (1 + ε )dxy −(1 + ε )xyz + xyz − bz2

 T  ≤

x y z

− ( 1 + ε )|d | 0 0

0 −ε lx 0

x y z

− ( 1 + ε )|d | 0 0

0 −ε lx∗ 0

 T  ≤

 T x y z



 

x Q5 y z

 

0 0 −b

x y z

 

0 0 −b

x y z

< 0.

Then, it has λM (Q5 )

V5 (X ) ≤ V5 (X0 )e λM (P5 )

(t−t0 )

.

(15)

Thus, (y, z) converge to (0, 0) or into the area of y < 0, z ≤ ly exponentially. Part 3. When x ≥ 0, y ≤ 0, z ≥ −hy (in OCDA¯ 1 O), choose a positive definite Lyapunov function as follows:

1 − ε¯ 1 − ε¯ 2 1 2 |d |x2 + y + z 2a 2 2 ⎡ ⎤ 1 − ε ¯ ⎡ ⎤T ⎡ ⎤ ⎡ ⎤T ⎡ ⎤ | d | 0 0 x ⎢ 2a x x ⎥ x ⎢ ⎥ 1 − ε ¯ ⎣ ⎦ ⎣ ⎦ ⎣ ⎦ ⎣ = y ⎢ 0 0 ⎥ y  y P6 y⎦. 2 ⎣ ⎦ z z z z 1 0 0 2

V6 =

(16)

Let ε¯ ∈ (0, 1 ) and x ≥ x∗ > 0, the derivative of V6 along system (1) is given by:

dV6 = −(1 − ε¯ )|d|x2 + (1 − ε¯ )|d|xy + (1 − ε¯ )dxy dt −(1 − ε¯ )xyz + xyz − bz2 ≤ −(1 − ε¯ )|d|x2 + ε¯ xyz − bz2

⎡ ⎤T ⎡ x

≤ ⎣y⎦ ⎣

0

⎤⎡ ⎤ x

0 ⎦⎣y⎦

0

−xε¯ h

0

0

−b

x

−(1 − ε¯ )|d|

0

0

0

−x∗ ε¯ h

0

0

0

−b

⎣y⎦ ⎣ z

 T 

0

z

⎡ ⎤T ⎡ ≤

−(1 − ε¯ )|d|

x y z

⎡ ⎤

z

⎤⎡ ⎤ x

⎦⎣y⎦ z

x

Q6 ⎣y⎦ < 0. z

Then, one has λM (Q6 )

V6 (X ) ≤ V6 (X0 )e λM (P6 )

(t−t0 )

.

(17)

Thus, (y, z) exponentially converges to (0, 0) or into y > 0, z ≤ −hy. But when (y, z) go into y > 0, z ≤ −hy or y < 0, z ≤ ly, by using Lyapunov function V4 , it can control (y, z) tends to (0, 0) exponentially again. Furthermore, by using Lemma 1, x also converges to 0 exponentially. Consequently, S0 is globally exponentially stable if d < 0.

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

215

Necessity. If S0 is globally exponentially stable, linearized matrix A0 should be a Hurwitz matrix, that is, max Reλ(A0 ) < 0. Since

max Reλ(A0 ) < 0 ⇔ max Reλ and

 λ + a   −d ⇒ λ1,2 =

−a

a

d

0



< 0,

  = λ2 + λa − ad = 0 λ

−a

√ a2 + 4ad . 2

−a ±

So, Re(λ1, 2 ) < 0⇔d < 0. The proof is complete.



Theorem 3. Sufficient and necessary condition for global asymptotic stability of S0 of system (1) is d = 0. Proof. Sufficiency. If d = 0, the second and third equations of system (1) become

⎧ dy ⎪ ⎨ = −xz, dt

(18)

⎪ ⎩ dz = xy − bz. dt

Here, we consider x as a parameter and construct a Lyapunov function as follows:

V7 =

1 2 1 2 y + z . 2 2

(19)

The derivative of V7 along system (18) is given by:

dV7 = −bz2 ≤ 0. dt

(20)

By part variable stability theorem, it is easy to obtain that the solution of (0, 0) for system (18) is stable in sense of Lyapunov, furthermore, z converges to zero asymptotically. dV Let dt7 |z(0 ) = 0, one has z = 0. By LaSalle’s invariance principle [30], it has

dy = −x × 0 = 0. dt Then, y(t) is a constant, let y = η and plug y = η into the first equation of system (1), it yields

dx = −ax + aη. dt Thus,

x(t ) = x0 e−a(t−t0 ) +



t

t0

by L’hospital’s rule,

lim x(t ) = lim

t→+∞

x0 +

e−a(t−τ ) aηdτ , t t0

ea(t−t0 ) aηdτ

ea(t−t0 )

t→+∞

= η.

Plugging x = y = η into the third equation of system (1), one obtains

z(t ) = z0 e−b(t−t0 ) +



t

t0

e−b(t−τ ) η2 dτ .

Similarly, it holds that

lim z(t ) =

t→+∞

η2 b

.

Since lim z(t ) = 0, so, η = 0. Plugging η = 0 into the first equation of system (1), one has t→+∞

lim x(t ) = η = 0.

t→+∞

Thus, S0 is globally asymptotically stable. Sufficiency is proved. Necessity. The necessity of asymptotic stability of system (1) on S0 is that max Reλ(A0 ) ≤ 0. From the proof of Theorem 2, one can see that max Reλ(A0 ) ≤ 0⇔d ≤ 0. As it is proved that the sufficient and necessary condition of global exponential

216

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

stability of system (1) on S0 is d < 0, then, the necessity of global asymptotic stability of system (1) on S0 is d = 0. The proof is complete.  Theorem 4. Sufficient and necessary condition for instability of system (1) on S0 is d > 0. Proof. Sufficiency. Three eigenvalues of A0 is

λ1 = −b, λ2 =

−a −

√ √ a2 + 4ad −a + a2 + 4ad , λ3 = . 2 2

Re(λ3 ) > 0 if d > 0. Thus, the linearized system of system (1) on S0 is unstable. According to the first order approximation theory, the corresponding nonlinear system (1) on S0 is unstable. Necessity. It is proved that sufficient and necessary condition of global exponential stability of system (1) on S0 is d < 0, sufficient and necessary condition of global asymptotic stability of system (1) on S0 is d = 0, then, the necessary condition of instability of system (1) on S0 is d > 0. The proof is complete.  Remark 2. Sufficient and necessary condition of global exponential stability, global asymptotic stability and instability for Yang–Chen system on S0 are studied from Theorem 2 to 4. There is nothing to do with the system when it is unstable. So, it is important to add a negative feedback controller to the system such that the closed loop system is globally stable. Furthermore, if the controller is simple and easily to be applied in reality, it is a good controller. For this purpose, we propose the following linear feedback controller to stabilize system (1) when it is unstable. Theorem 5. When d > 0, a negative feedback controller −ky is added into the second equation of system (1). If k = d, then the feedback controlled system (1) is globally asymptotically stable on S0 ; if k = d + ε with ε > 0, then the feedback controlled system (1) is globally exponentially stable on S0 . Proof. Add the feedback controller u = −ky to the second equation of system (1), it yields

⎧ dx ⎪ = a ( y − x ), ⎪ ⎪ dt ⎪ ⎪ ⎨ dy

⎪ dt ⎪ ⎪ ⎪ ⎪ ⎩ dz dt

= dx − xz − ky,

(21)

= xy − bz.

When k = d + ε , construct a radially unbounded and positive definite Lyapunov function

⎡ ⎤T ⎡ ⎤

V8 =

x x x y z + +  ⎣y⎦ P8 ⎣y⎦. 2a 2d 2d z z 2

2

2

(22)

Take the derivative of V8 along (21), it holds that

d+ε 2 dV8 = −x2 + xy + xy − y − dt d ⎡ 1  T −1 x ⎢ ε −1 − = y ⎢1 d ⎣ z 0 0

⎡ ⎤T x



⎡ ⎤

b 2 z d ⎤ 0 ⎡

⎤ x ⎥ 0 ⎥⎣y⎦ ⎦

−

b d

z

x

⎣y⎦ Q8 ⎣y⎦ z

z

⎡ ⎤T ⎡ ⎤ x

x

z

z

λM (Q8 ) ≤ λM (Q8 )⎣y⎦ ⎣y⎦ ≤ V . λm (P8 ) 8 Then, it has λM (Q8 )

V8 (X ) ≤ V8 (0 )e λm (P8 )

(t−t0 )

.

(23)

Hence, S0 is the global exponential stable equilibrium of system (21). If k = d, one has

b dV8 |(23) = −(x − y )2 − z2 ≤ 0. dt d

(24)

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

217

So, by the part variable stability theorem, x − y and z converge to zero globally asymptotically. It is easy to see from the third equation of system (21) that x and y also converge to zero asymptotically.  Corollary 1. When d > 0, a negative feedback controller −kx is added into the second equation of system (1). If k = d, the feedback controlled system (1) is globally asymptotically stable on S0 ; if k = d + ε with ε > 0, the feedback controlled system (1) is globally exponentially stable on S0 . Proof. With the feedback controller, system (1) becomes

⎧ dx ⎪ = a ( y − x ), ⎪ ⎪ dt ⎪ ⎪ ⎨ dy = (d − k )x − xz, ⎪ dt ⎪ ⎪ ⎪ ⎪ ⎩ dz = −bz + xy.

(25)

dt

According to Theorem 2 and 3, it is easy to obtain the result.



5. Sufficient and necessary condition and application for Lyapunov stability of equilibrium points of S+ , S− As it is discussed in Section 4, system (1) is globally stable if d ≤ 0 and S0 is the equilibrium point. If d > 0, S0 is not a stable equilibrium anymore. Since there exist other two equilibrium points S+ , S− , this section focus on the study of local stability properties of these two equilibria. Theorem 6. The characteristic polynomial of A+ (A− ) is given by:

⎡

λ+a

−a

0 √ − bd

√ − bd

f (λ ) = det(λ − A+ ) = det⎣

λ

⎤

0 √ bd ⎦

λ+b = λ + (a + b)λ + (ab + bd )λ + 2abd = det(λ − A− ). 3

2

And d > 0. Then, the following conditions are equivalent. Case 1. f(λ) is a Hurwitz polynomial, +b ) ⇔ a(aa−b > d, ⇔ the linearized system of system (1) on S+ (S− ) is locally exponentially stable, that is, system (1) is locally exponentially stable on S+ (S− ). Case 2. The coefficient matrix A+ (A− ) has at least one eigenvalue with positive real part, ⇔ the linearized system of system (1) on S+ (S− ) is unstable, +b ) ⇔ a(aa−b < d. Case 3. The linearized systems of system (1) on S+ (S− ) is locally stable while non-asymptotically stable, +b ) ⇔ a(aa−b = d. Proof. Case 1. The former two equivalent conditions are in fact the special case of n = 3 in Hurwitz theorem (see Theorem 2.6.2 in [21]). It remains to prove the third equivalent condition only. ⇒ It is an immediate corollary of Theorem 9.1 in [22]. ⇐ Assuming (x∗ , y∗ , z∗ ) is one of the two equilibrium points S+ , S− . Let x¯ = x − x∗ , y¯ = y − y∗ , z¯ = z − z∗ , one has

⎧ ¯ dx ⎪ ⎪ dt = a(y¯ − x¯ ), ⎪ ⎪ ⎪ ⎨ dy¯

⎪ dt ⎪ ⎪ ⎪ ⎪ ⎩ dz¯ dt

= −x¯z¯ − x∗ z¯,

(26)

= x¯y¯ + y∗ x¯ + x∗ y¯ − bz¯.

The corresponding linearized system is

⎧ ¯ dx ⎪ ⎪ dt = a(y¯ − x¯ ), ⎪ ⎪ ⎪ ⎨ dy¯

⎪ dt ⎪ ⎪ ⎪ ⎪ ⎩ dz¯ dt

= −x∗ z¯,

= y∗ x¯ + x∗ y¯ − bz¯.

(27)

218

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

Since system (26) is locally exponentially stable on S+ (S− ), in a neighborhood of S+ (S− ), there exist two positive constants M > 0, λ > 0 such that

|x¯| ≤ Me−λ(t−t0 ) , |y¯| ≤ Me−λ(t−t0 ) , |z¯| ≤ Me−λ(t−t0 ) . Let 0 < ε < λ2 , and take the following transformation:

ξ = eεt x¯, η = eεt y¯, γ = eεt z¯. System (26) is rearranged as follows:

⎧ ⎨ξ˙ = (−a + ε )ξ + aη, η˙ = −ξ γ e−εt − εη − x∗ γ , ⎩ γ˙ = ξ ηe−εt + y∗ ξ + x∗ η + (−b + ε )γ .

(28)

It is easy to obtain that

λ − (t − t0 ) |ξ (t )| ≤ |eεt ||x¯ (t )| ≤ Meεt e−λ(t−t0 ) ≤ Me 2 , λ − (t − t0 ) |η (t )| ≤ |eεt ||y¯ (t )| ≤ Meεt e−λ(t−t0 ) ≤ Me 2 , λ − (t − t0 ) |γ (t )| ≤ |eεt ||z¯ (t )| ≤ Meεt e−λ(t−t0 ) ≤ Me 2 . So, the equilibrium point (0, 0, 0) of system (28) is also exponentially stable. For system (28), the corresponding linearized system on (0, 0, 0) is of the following form:

⎧ ⎨ξ˙ = (−a + ε )ξ + aη, η˙ = −εη − x∗ γ , ⎩ γ˙ = y∗ ξ + x∗ η + (−b + ε )γ .

(29)

The coefficient matrix of system (27) is as follows:

⎡

−a

a

A=⎣ 0

0

∗

∗

y

x

0

⎤

−x∗ ⎦, −b

and the coefficient matrix of system (29) is given by:



A¯ =

−a + ε 0 y∗

a

ε

x∗



0 −x∗ . −b + ε

From (26) and (28), it is derived that the solution is exponentially stable. Moreover, max Reλ(A) < 0 and max Reλ(A¯ ) < 0. Otherwise, there will come forth a contradiction with the first-order approximation instability theory. Since max Reλ(A ) = max Reλ(A¯ ) − ε < −ε , the linearized system corresponding to (27) is also exponentially stable. Case 2. For the cubic polynomial f(λ) with real coefficients, there is no doubt that there exists at least one real root. Furthermore, this real root is bound to be negative in the event that all the coefficients are positive. Assuming λ1 = α1 < 0, λ2,3 = α ± iβ , one has

λ1 + λ2 + λ3 = α1 + 2α = −(a + b), λ1 λ2 + λ1 λ3 + λ2 λ3 = 2α1 α + α 2 + β 2 = ab + bd, λ1 λ2 λ3 = α1 (α 2 + β 2 ) = −2abd, and

(a + b)(ab + bd ) − 2abd = (−α1 − 2α )(2α1 α + α 2 + β 2 ) + α1 (α 2 + β 2 ) = −2α (α1 + α )2 − 2αβ 2 < 0 ⇔ α > 0. As a result, (a + b)(ab + bd ) − 2abd < 0 ⇔ the coefficient matrix of the linear system corresponding to S+ (S− ) has an eigenvalue with a positive real part ⇔ the linear system is unstable of the positive Lyapunov exponential type ⇔ the nonlinear system corresponding to S+ (S− ) is unstable of the positive Lyapunov exponential type. Case 3. It follows from the proof of Case 2 that (a + b)(ab + bd ) − 2abd = 0 ⇔ α = 0 ⇔ there is a negative real root and two conjugate pure imaginary roots of f (λ ) = 0. So, the corresponding linearized system on S+ (S− ) is Lyapunov stable while non-asymptotically stable. The stability problem of the nonlinear system corresponding to S+ (S− ) belongs to the second-critical state. That will be our future topic. 

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

219

Fig. 4. The trajectory X(t) of Yang–Chen system with d < 0.

Fig. 5. The trajectory X(t) of Yang–Chen system with d = 0.

6. Discussing the question about bifurcation The bifurcation is one of the core problems in the qualitative analysis of differential dynamical systems, especially for chaos control and anti-control. If a least negative feedback controller is added near the bifurcation, the system structure can be changed. We discuss the bifurcation of system (1) in the following. (I) The branch solution of equilibrium S0 for system (1) is d = 0 on condition that a, b are positive constants, since system (1) is globally exponentially stable if d < 0, it is globally asymptotically stable if d = 0, and it is unstable if d > 0. (II) As for S+ (S− ), not only d but also a, b may result in bifurcations, and a, b, d are distributed on a hypersurface in a canonical parameter space. (a) ∀a, b > 0, a > b, a small change of d may lead to the following results. +b ) +b ) If d < a(aa−b , S+ (S− ) is local exponentially stable; if d = a(aa−b , S+ (S− ) is local stable while non-asymptotically stable; +b ) if d > a(aa−b , S+ (S− ) is unstable.

+b ) Thus, d = a(aa−b is a branch solution of equilibrium S+ (S− ) for system (1) if ∀a, b > 0 and a > b. √ d−b± (d−b)2 −4ad (b) Similarly, a = is a branch solution of equilibrium S+ (S− ) for system (1) if ∀a, b > 0, d > 0, d > b, (d − 2 b)2 ≥ 4ad. ) (c) b = a(ad−a is also a branch solution of equilibrium S+ (S− ) for system (1) if ∀a, b > 0, a > b, d > a. +d

7. Simulation In this section, some simulation results will be proposed to illustrate the theoretical results obtained in the previous sections. A fourth order Runge–Kutta method is used to obtain the simulation results with MATLAB software. By using different initial states and parameters, the figures are obtained as follows. Let a = 10, b = 1, if d = −1 and initial state is [x0 , y0 , z0 ] = [1, −1, −2], Fig. 4 shows that the equilibrium S0 of system (1) is globally exponentially stable. If d = 0 with initial state [x0 , y0 , z0 ] = [1, −1, −2], Fig. 5 shows the good result of global asymptotic stability of system (1). If d > 0, +b ) even d is very small, i.e., d = 0.01, equilibrium S0 is not global stable any more, since a(aa−b = 11 9 > 0.01 = d, but equilibrium

220

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

Fig. 6. The local exponential stability of of Yang–Chen system at S+ with d > 0.

Fig. 7. The local exponential stability of Yang–Chen system at S− with d > 0.

Fig. 8. The local stability of Yang–Chen system at S+ with d > 0.

S+ (S− ) is locally exponentially stable, Figs. 6 and 7 show the good simulation results. If let d = but non-asymptotically stable, Figs. 8 and 9 show the simulation results.

11 9 ,

S+ (S− ) is locally stable,

8. Conclusion This paper has studied the Lagrange stability and Lyapunov stability of Yang–Chen system. By introducing different radial unbounded Lyapunov functions in different regions, a global exponential attractive set of Yang–Chen chaotic system has been constructed with geometrical and mathematical methods. Then, simple algebraic necessary and sufficient conditions of global exponential stability, global asymptotic stability, and global instability of equilibrium S0 (0, 0, 0) were proposed. And we also obtained the relevant expression of the corresponding parameters for local exponential stability, local asymptotic stability, local instability of equilibria S± . Furthermore, some branch expressions were given for the parameters of the system.

X. Liao et al. / Applied Mathematics and Computation 309 (2017) 205–221

221

Fig. 9. The local stability of Yang–Chen system at S− with d > 0.

Lagrange stability and Lyapunov stability conditions presented here are very useful for studying how chaos produced, chaos control, chaos synchronization, and other applications. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36]

E.N. Lorenz, Deterministic non-periodic flow, J. Atoms. Sci. 20 (1963) 130–141. E.N. Lorenz, The Essence of Chaos, University of Washington Press, Seattle, 1993. C. Sparrow, The Lorenz Equations: Bifucation, Chaos and Strange Attactors, Springer, New York, 1982. O. Rössler, An equation for continuous chaos, Phys. Lett. A 57 (1976) 397. L.O. Chua, The genesis of Chua’s circuit, Arch. Elektron. Ubertragungstechnik 46 (4) (1992) 250–257. G.R. Chen, X.N. Dong, From Chaos to Order, Word Scientific Publishing Co., Singapore, 1998. I. Stewart, The Lorenz attractor exists, Nature 400 (2002) 948–949. J.H. Lü, T.S. Zhou, G.R. Chen, et al., The compound structure of Chen’s attractor, Int. J. Bifurc. Chaos 12 (2002) 855–858. G.R. Chen, J.H. Lü, Control and synchroniztion of Lorenz familics, Dyn. Anal. (2003a). Chinese Press G.R. Chen, J.H. Lü, Chen attractor exists, Int. J. Bifurc. Chaos 14 (2003b) 3167–3177. Q. Ma, S. Xu, F.L. Lewis, B. Zhang, Y. Zou, Cooperative output regulation of singular heterogeneous multiagent systems, IEEE Trans. Cybern. 46 (6) (2016) 1471–1475. Q. Ma, Cooperative control of multi-agent systems with unknown control directions, Appl. Math. Comput. 292 (2017) 240–252. Q.G. Yang, G.R. Chen, A chaotic system with one saddle and two stable node-foci, Int. J. Bifurc. Chaos 18 (5) (2008) 1393. G. Leonov, A. Bunin, N. Koksch, A tractor localization of the Lorenz system, ZAMM 67 (1987) 649–656. G. Leonov, On estimates of attractors of Lorenz system, Vestn Leningr. Unuversiten Mat. 21 (1988) 32–37. G. Leonov, Bound for attrators of the existence of homoclinic orbits in the Lorenz system, J. Appl. Math. 65 (2001) 19032. Q. Ma, S. Xu, Y. Zou, G. Shi, Synchronization of stochastic chaotic neural networks with reaction–diffusion terms, Nonlinear Dyn. 67 (3) (2012) 2183–2196. X.X. Liao, B.J. Xu, P. Yu, G.R. Chen, Constructive proof of globally exponentially attractive and positively invariant set of the chaotic Chen’s system, Sci. China Inf. Sci. 45 (1) (2015) 129–144. S.C. Jeong, D.H. Ji, J.H. Park, S.C. Won, Adaptive synchronization for uncertain chaotic neural networks with mixed time delays using fuzzy disturbance observer, Appl. Math. Comput. 219 (11) (2013) 5984–5995. S.J. Choi, S.M. Lee, S.C. Won, J.H. Park, Improved delay-dependent stability criteria for uncertain lure systems with sector and slope restricted nonlinearities and time-varying delays, Appl. Math. Comput. 208 (2) (2009) 520–530. S.M. Lee, J.H. Park, Robust stabilization of discrete-time nonlinear lure systems with sector and slope restricted nonlinearities, Appl. Math. Comput. 200 (15) (2008) 429–436. X.X. Liao, Y.L. Fu, S.L. Xie, et al., Globally exponentially attractive sets of the family of Lorenz systems, Sci. China Ser. F-Inf. Sci. 51 (2008) 283–292. X.X. Liao, On new results of global attractive set and positively invariant sets for the Lorenz chaotic system and the application in chaos control and synchronization, Sci. Sin. Technol. 34 (2004) 1404–1419. X.X. Liao, H.G. Luo, Y.L. Fu, et al., Analysis on the global exponential set and positive invariant set of the Lorenz family, Sci. Sin. Technol. 37 (2007) 757–769. Q. Luo, X.X. Liao, Simple algebraic necessary and suffcient conditions of Lyapunov stablility for the Lorenz chaotic system and its application, Sci. Sin. Technol. 40 (2010) 1086–1095. G.P. Zhou, J.H. Huang, X.X. Liao, S.J. Cheng, Stability analysis and control of a new smooth Chua’s system, Abstr. Appl. Anal. 2013 (620286) (2013) 1–10. G.P. Zhou, D. Liu, J.H. Huang, X.X. Liao, Global synchronization of a new Chua’s system, Int. J. Bifurc. Chaos 25 (4) (2015) 1550053(1–13). X.X. Liao, B.J. Xu, P. Yu, et al., Constructive proof of globally exponentially attractive and positively invariant set of the chaotic Chen’s system, Sci. Sin. Technol. 34 (2004) 1404–1419. L. Solomon, Stability of Nonlinear Control Systems, Academic Press, New York, 1965. J.P. Lasalle, The Stability of Dynamical Systems, SIAM Philadelphia, 1976. E. Tlelo-Cuautle, L.G. de la Fraga, J.J.R. Magdaleno, Engineering Applications of FPGAs: Chaotic Systems, Artificial Neural Networks, Random Number Generators, and Secure Communication Systems, Springer, 2016, doi:10.1007/978- 3- 319- 34115- 6. A.D. Pano-Azucena, J.J. Rangel-Magdaleno, et al., Arduino-based chaotic secure communication system using multi-directional multiscroll chaotic oscillators, Nonlinear Dyn. (2016), doi:10.1007/s11071- 016- 3184- 4. E. Tlelo-Cuautle, A.J. Quintas-Valles, L.G. de la Fraga, et al., VHDL descriptions for the FPGA implementation of PWL-function-based multi-scroll chaotic oscillators, PLoS ONE 11 (12) (2016a) 1–32. E. Tlelo-Cuautle, A.D. Pano-Azucena, J.J. Rangel-Magdaleno, et al., Generating a 50-scroll chaotic attractor at 66 MHz by using FPGAs, Nonlinear Dyn. 85 (4) (2016b) 2143–2157. E. Tlelo-Cuautle, V.H. Carbajal-Gomez, P.J. Obeso-Rodelo, et al., FPGA realization of a chaotic communication system applied to image processing, Nonlinear Dyn. 82 (4) (2015a) 1879–1892. E. Tlelo-Cuautle, J.J. Rangel-Magdaleno, A.D. Pano-Azucena, et al., FPGA realization of multi-scroll chaotic oscillators, Commun. Nonlinear Sci. Numer. Simul. 27 (1–3) (2015b) 66–80.