Dynamical analysis of a new autonomous 3-D chaotic system only with stable equilibria

Nonlinear Analysis: Real World Applications 12 (2011) 106–118

Contents lists available at ScienceDirect

Nonlinear Analysis: Real World Applications journal homepage: www.elsevier.com/locate/nonrwa

Dynamical analysis of a new autonomous 3-D chaotic system only with stable equilibria Zhouchao Wei ∗ , Qigui Yang 1 School of Mathematical Sciences, South China University of Technology, Guangzhou, 510640, PR China

article

info

Article history: Received 17 January 2010 Accepted 18 May 2010 Keywords: Chaotic attractors Degenerate heteroclinic cycles Sil’nikov theorem Lyapunov exponent Poincaré map

abstract This paper presents a new 3-D autonomous chaotic system, which is topologically nonequivalent to the original Lorenz and all Lorenz-like systems. Of particular interest is that the chaotic system can generate double-scroll chaotic attractors in a very wide parameter domain with only two stable equilibria. The existence of singularly degenerate heteroclinic cycles for a suitable choice of the parameters is investigated. Periodic solutions and chaotic attractors can be found when these cycles disappear. Finally, the complicated dynamics are studied by virtue of theoretical analysis, numerical simulation and Lyapunov exponents spectrum. The obtained results clearly show that the chaotic system deserves further detailed investigation. © 2010 Elsevier Ltd. All rights reserved.

1. Introduction In 1963, Lorenz found the first chaotic attractor in a simple mathematical model of a weather system which was made up of three linked nonlinear differential equations [1]. In 1976, Rossler conducted important work that rekindled the interest in low-dimensional dissipative dynamical systems [2]. Sprott embarked upon an extensive search for autonomous threedimensional chaotic systems with fewer than seven terms in the right-hand side of the model equations. He considered general three-dimensional ordinary differential equations with quadratic nonlinearities and found by computer simulation 19 simple 3-D quadratic autonomous chaotic systems with none, one equilibrium or two equilibria [3–5]. Later, many Lorenz-like or Lorenz-based chaotic systems were proposed and investigated. Some classical 3-D autonomous chaotic systems have three particular fixed points: one saddle and two unstable saddle–foci [1,6–8]. The other 3-D chaotic systems have two unstable saddle–foci [2,9,10]. In 2008, Yang and Chen found another 3-D chaotic system with three fixed points: one saddle and two stable equilibria [11]. Moreover, many theoretical analysis and numerical simulation results about these systems are obtained [12–23]. It should be noted that one commonly used analytic criterion for generating and proving chaos in autonomous systems is based on the fundamental work of Sil’nikov [24,25] and its subsequent embellishment and slight extension [26]. All of these outstanding results require the system to have at least an unstable equilibrium. However, Sil’nikov criteria is sufficient but certainly not necessary for emergence of chaos. Therefore, it is interesting to ask whether or not there are 3-D autonomous chaotic systems only with stable equilibria and at most six terms including only one or two nonlinear terms such that the topological structure of the chaotic systems and the well-known 3-D autonomous chaotic systems are different? In this paper, this question has been answered. We introduce a new 3-D chaotic system with six terms including only one nonlinear term in the form of exponential function. It is surprising to find when all of equilibria are stable, the system generates a double-scroll chaotic attractor, not satisfying the criteria in the traditional Sil’nikov homoclinic/homoclinic theorem.

∗

Corresponding author. Tel.: +86 20 85294419. E-mail addresses: [email protected] (Z. Wei), [email protected] (Q. Yang).

1 Tel.: +86 20 87110448. 1468-1218/$ – see front matter © 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.nonrwa.2010.05.038

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

107

Table 1 Numerical analysis for several typical chaotic attractors of system (1) with parameters (a, m, n) = (0.8696, 0.756144, 10.5) and varying b at initial values (−0.08, −1.02, −0.49). Parameter b

Equilibria

Eigenvalues

LE

Fractional dimensions DL = 2.259315

b = 0.5

(±1.53342, ±1.53342, 0.66125)

λ1 = −1.69015 λ2,3 = 0.160275 ± 4.38005i

λL1 = 0.4804 λL2 = −0.0009 λL3 = −1.8491

b = 0.8696

(±1.53342, ±1.53342, 1.15005)

λ1 = −1.7392 λ2,3 = ±4.32074i

λL1 = 0.3606 λL2 = −0.0005 λL3 = −2.0993

DL = 2.171533

b = 2.1

(±1.53342, ±1.53342, 2.77725)

λ1 = −1.94684 λ2,3 = −0.511378 ± 4.05168i

λL1 = 0.2456 λL2 = −0.0006 λL3 = −3.2146

DL = 2.076215

The complicated dynamics are analyzed by theoretical analysis, numerical simulation and Lyapunov exponents spectrum. The evolution processes of this system and dynamics behaviors will be presented when parameters and initial values vary. The system can produce chaotic attractors coexistence with three types of equilibria: saddle–foci, non-hyperbolic and stable. By the calculation of the first Lyapunov coefficient, we determine the stability of the two non-hyperbolic equilibria, and there are limit cycles bifurcating from them. In particular, when the parameter values are fixed and small changes in initial values, the chaotic attractor can also be obtained in the asymptotically stable region of all of the equilibria. It is immediately clear that the system will be topologically non-equivalent to the original Lorenz and all Lorenz-like systems. In addition, We present the existence of singularly degenerate heteroclinic cycles [16,17] and numerically find periodic solutions and chaotic attractors when these cycles disappear. It is organized as follows. In Section 2, we present the new system and report the existence of chaos with different types of equilibria and initial conditions. In Section 3, we analyze stability of the equilibria, and a brief review of the methods used to study codimension one Hopf bifurcation are presented. Simulation results are then presented in Section 4, which visualize and illustrate the very rich complex dynamical behaviors with some interesting characteristics. Also, we numerically find periodic solutions and chaotic attractors when the singularly degenerate heteroclinic cycles disappear. Finally, in Section 5, we make some concluding remarks. 2. The chaotic system In this section, we present a new 3-D autonomous ODE system: x˙ = a(y − x) y˙ = −by + mxz z˙ = n − exy ,



(1)

√

√

where and n are real parameters. If n > 1 and m ̸= 0, it possesses two equilibria E1 ( ln n, ln n, b/m) and √ a, b , m √ E2 (− ln n, − ln n, b/m). The system has the simple form, but can display complicated and unusual dynamical behaviors. We will show that, by traditional several ways of describing chaotic attractor, the system is indeed chaotic. When parameters a = 0.8696, m = 0.756144, n = 10.5 and b varies, we analyzes the evolution processes of this system and dynamics behaviors with initial values (−0.08, −1.02, −0.49). For b = 0.5, it displays a double-scroll chaotic attractor with novel Lorenzunlike shaped, as shown in Fig. 1(a), (d). Moreover, We observe that system (1) has a chaotic attractor for parameter b = 0.8696 (see Fig. 1(b), (e)), coexisting with two non-hyperbolic equilibria whose characteristic values both are: λ1 = −1.7392, λ2,3 = ±4.32074i. Therefore, equilibria E1,2 have local one-dimensional stable manifolds tangent to the eigenvector associated to the eigenvalue λ1 and local two-dimensional center manifolds tangent to the plane generated by the eigenvectors associated to the complex eigenvalues λ2,3 . In fact, the equilibria E1,2 are both unstable by calculating the first Lyapunov coefficient (Section 3). In particular, for parameter b = 2.1, another type of chaotic attractor can be observed (see Fig. 1(c), (e)). The Poincaré map also shows that this system is chaotic. It is noted that, the chaotic attractor is a bit different from the chaotic attractor plotted in Fig. 1(a), because the equilibria E1,2 both are stable, whose characteristic values both are λ1 = −1.94684, λ2,3 = −0.511378 ± 4.05168i. Therefore, system (1) has neither homoclinic orbits nor heteroclinic orbits joining E1 and E2 . Of particular interest is that system (1) can display a chaotic attractor, not satisfying the Sil’nikov theorem (see Fig. 1(c), (e)). In addition, the shapes of the two chaotic attractors are different from that of the Lorenz system or any existing systems. Furthermore, we examine three Lyapunov exponents of the system (1), which are denoted by Li for i = 1, 2, 3 and ordered as L1 > L2 > L3 . Note that the system (1) has a chaotic attractor if L1 > 0, L2 = 0, L3 < 0. To clearly show the above three types of chaotic attractors mentioned previously, we list equilibria, eigenvalues, Lyapunov exponents and fractional dimensions for the typical parameter values in Table 1. It is noted that, for b = 2.1, all equilibria of the system are stable, chaotic attractor has been obtained with initial values (−0.08, −1.02, −0.49).

108

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

a

b

20

40 20

0 z

z

0 –20

–20

–40 –40

–60 10

5

2 0

2 0

0

y –5

y

x

–2

–10

c

d

x

–2

6

40

4

20

2

0 y

z

0

0

–20 –2 –40 5 0

y –5

e

–6

x

–2

f

6

–2

–1

0 x

1

2

–6 –2

–1

0 x

1

2

6

4

4

2

2

0

0

y

y

–4

2 0

–2

–2

–4

–4

–6

–2

–1

0 x

1

2

Fig. 1. Chaotic attractors and Poincaré maps of system (1) with (a, m, n) = (0.8696, 0.756144, 10.5) at the initial values (−0.08, −1.02, −0.49): (a) b = 0.5; (b) b = 0.8696; (c) b = 2.1; (d) Poincaré map of (a) on plane z = −5; (e) Poincaré map of (b) on plane z = −5; (f) Poincaré map of (c) on plane z = −5.

Now, we investigate the influence of initial condition on the dynamics of system with the fixed parameters values. In particular, when we fix a = 0.8696, b = 2.1, m = 0.756144, n = 47 and change initial values slightly, dynamical behaviors of the system may produce large variations in the long term. Besides the two stable equilibrium points, chaotic attractors of system (1) also are obtained, which imply that chaos coexists with the two stable fixed points: (a) A chaotic attractor with initial values (−0.08, −3.65, −0.41) is shown as Fig. 2(a). The Lyapunov exponents of the system are found to be L1 = 0.6706, L2 = −0.0003, and L3 = −3.6399; (b) For initial values (−0.09, −3.65, 0.41), trajectories converge to stable equilibrium E1 (see Fig. 2(b)), The Lyapunov exponents of the system are found to be L1 = −0.6002, L2 = −0.6015 and L3 = −1.7679; (c) A chaotic attractor with with initial values (−0.09, −3.62, −0.41) is shown as Fig. 2(c). The Lyapunov exponents of the system are found to be L1 = 0.6784, L2 = −0.0002, and L3 = −3.6478. A small change in the initial condition of the system causes wide difference of trajectories. For different initial conditions, trajectories converge to different attractors. The above investigations verify that when all of equilibria of system (1) are stable, trajectory can converges to two types of attractors (chaos or stable equilibrium). z-coordinate (t , z (t )) of the three attractors show the sensitive dependence on initial values (see Fig. 2(d)).

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

b

c

2

10

5

0

5

0

2

0

y

y

10

–5

y

a

–4

–10 100

–5

–6 20 50

2

–10 100 10

0

1

0 z

–50 –100 –2

d

0 –1

z

x

109

0 –10 –20 –2

50

2 1

0

–1

z

x

–50 –100 –2

0 –1

x

80

50

z(t)

20

–10

–40

–70

–100

0

20

40

60 t

80

100

120

Fig. 2. Phase portraits of system (1) with a = 0.8696, b = 2.1, m = 0.756144 and n = 47 for three sets of initial values: (a) initial values (−0.08, −3.65, −0.41), (b) initial values (−0.09, −3.65, −0.41), (c) initial values (−0.09, −3.62, −0.41), (d) z-coordinate (t , z (t )) show the sensitive dependence on these initial values.

3. Some basic properties and bifurcation analysis of the new system (1) 3.1. Symmetry and invariance Firstly, it is easy to see the invariance of system under the coordinate transformation (x, y, z ) → (−x, −y, z ), i.e., the system has rotation symmetry around the z-axis. 3.2. Dissipativity For system (1), we can obtain the divergence

∇ ·V =

∂ x˙ ∂ y˙ ∂ z˙ + + = −(a + b). ∂x ∂y ∂z

This means system (1) is dissipative system when a + b > 0. Previous numerical simulations seem to suggest that solutions of the system are bounded. However, this is not true. Taking an initial condition (0, 0, k) and the parameter n ̸= 1, we have the solution of this system (0, 0, t + k), which is not bounded as the Rikitake system [27–29]. Thus it can not be in the basin of attraction of any chaotic attractor. 3.3. Equilibria and stability Note that n < 1 or n ̸= 1, m = 0, there is no equilibria. When n = √ 1, all the ), z ∈ R, are √ equilibria in the form √ (0, 0, z√ non-isolated. When n > 1 and m ̸= 0, system (1) has two equilibria: E1 ( ln n, ln n, b/m) and E2 (− ln n, − ln n, b/m). In following, we consider stability of the equilibria if n > 1 and m ̸= 0.

110

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

Due to symmetry, we only consider E1 and the Jacobian matrix of system (1) is given as follows:

 J (E1 ) = 

−a

a

b √ −n ln n

− √b −n ln n

√0



m ln n .

(2)

0

Obviously, the characteristic equation about the equilibria E1 is: det (λI − J (E1 )) = λ3 + (a + b)λ2 + mn ln nλ + 2amn ln n = 0.

(3)

According to the Routh–Hurwitz criterion, the real parts of all the roots λ are negative if and only if

∆1 = a + b > 0 ,    a+b 1  ∆2 =  = (b − a)mn ln n > 0, 2mna ln n mn ln n ∆3 = 2amn ln n ∆2 > 0. From these inequalities, there are m > 0,

n > 1,

b > |a|.

(4)

Therefore, the two equilibria E1 and E2 are both asymptotically stable when the above conditions are met. If b = a, the Jacobian matrix has a pair of purely imaginary eigenvalues and a non-zero real eigenvalue. Therefore, E1 and E2 are both not hyperbolic but weak repelling focus (see following analysis). 3.4. Bifurcation analysis in system (1) Now we introduce the projection method [18,19] for the calculation of the first Lyapunov coefficient associated to the Hopf bifurcation, denoted by l1 . Consider the differential equation X˙ = f (X , µ),

(5)

where X ∈ R and µ ∈ R are respectively vectors representing phase variables and control parameters. Assume that f is a class of C ∞ in R3 × R4 . Suppose that (5) has an equilibrium point X = X0 at µ = µ0 , and denoting the variable X − X0 also by X , write 3

4

F (X ) = f (X , µ0 ),

(6)

as 1 1 B(X , X ) + C (X , X , X ) + O(‖ X ‖4 ), 2 6 where A = fx (0, µ0 ) and, for i = 1, 2, 3, F (X ) = AX +

 3 − ∂ 2 Fi (ξ )  B(X , Y ) =  ∂ξj ∂ξk  j,k=1

 ∂ 3 Fi (ξ )  C (X , Y , Z ) =  ∂ξj ∂ξk ∂ξl  j,k,l=1 3 −

Xj Yk , ξ =0

Xj Yk Zl . ξ =0

Suppose that A has a pair of complex eigenvalues on the imaginary axis: λ2, 3 = ±iw0 (w0 > 0), and these eigenvalues are the only eigenvalues with Re λ = 0. Let T c be the generalized eigenspace of A corresponding to λ2, 3 . Let p, q ∈ C 3 be vectors such that Aq = iw0 q,

AT p = −iw0 p,

⟨p, q⟩ = 1,

where A is the transposed of the matrix A. Any vector y ∈ T c can be represented as y = w q + w¯ ¯ q, where w = ⟨p, y⟩ ∈ C . The two-dimensional center manifold associated to the eigenvalues λ2, 3 can be parameterized by w and w ¯ , by means of an immersion of the form X = H (w, w) ¯ , where H : C 2 → R3 has a Taylor expansion of the form T

H (w, w) ¯ = wq + w¯ ¯q+

−

1

j!k! 2≤j+k≤3

hjk w j w ¯ k + O(|w|4 ),

with hjk ∈ C 3 and hjk = h¯kj . Substituting this expression into (6) we obtain the following differential equation Hw w ′ + Hw w ¯ ′ = F (H (w, w)), ¯ where F is given by (6). The complex vectors hjk are obtained by solving the system of linear equations defined by the coefficients of (6). Taking into account the coefficients of F , system (6), on the chart w for a central manifold, can be written as follows 1 w ˙ = iw0 w + G21 w|w|2 + O(|w|4 ), 2

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

111

with G21 ∈ C . The first Lyapunov coefficient can be written as l1 =

1 2

Re G21 ,

where G21 = ⟨p, C (q, q, q¯ ) + B(¯q, h20 ) + 2B(q, h11 )⟩. A Hopf bifurcation point (X0 , µ0 ) is an equilibrium point of (5) where the Jacobian matrix A only has a pair of purely imaginary eigenvalues λ2, 3 = ±iw0 (w0 > 0), and the other eigenvalue with non-zero real part. At a Hopf point a twodimensional center manifold is well defined, it is invariant under the flow generated by (5) and can be continued with arbitrary high class of differentiability to nearby parameter values. A Hopf point is called transversal if the parameter-dependent complex eigenvalues cross the imaginary axis with nonzero derivative. In a neighborhood of a transversal Hopf point with l1 ̸= 0 the dynamic behavior of the system (5), reduced to the family of parameter-dependent continuations of the center manifold, is orbitally topologically equivalent to the following complex normal form

w ′ = (η + iw)w + l1 w|w|2 , where w ∈ C , η, w and l1 are real functions having derivatives of arbitrary higher order, which are continuations of 0, w0 and the first Lyapunov coefficient at the Hopf point [18]. As l1 < 0 (l1 > 0) one family of stable (unstable) periodic orbits can be found on this family of manifolds, shrinking to an equilibrium point at the Hopf point. In the rest of this section we employs the three-dimensional Hopf bifurcation theory and applies symbolic computations to perform the analysis of parametric variations with respect to dynamical bifurcations. Because the system is invariant under the transformation (x, y, z ) → (−x, −y, z ), we only consider the bifurcation of the system (1) at E1 . Under the following linear transformation

 √ x1 = x − √ln n, y1 = y − ln n, z1 = z − b/m, which transforms the equilibrium E1 to the origin O(0, 0, 0), system (1) becomes

 x˙ 1 = a(y1 − x1 ), √ y˙ 1 = bx1 − by1 +√m ln n z1 + mx1 z1 ,  z˙1 = n − nex1 y1 + ln n (x1 +y1 ) .

(7)

The Jacobian matrix A of system (7) and its corresponding characteristic equation are same to (2) and (3), respectively. When (a, b, m, n) ∈ {(a, b, m, n)| √a > 0, b > 0, m > 0, n > 1, a = a0 = b}, (3) possesses a negative real root −2b and conjugate purely imaginary roots ± mn ln n i. Under this condition, the transversality condition

λ′ (a = a0 ) =

3mn ln n 9mn ln n + 16b2

>0

(8)

is also satisfied. Accordingly, Hopf bifurcation at E1 occurs. The stability of E1 depends on the value of the first Lyapunov coefficient l1 . We have the following theorem. Theorem 3.1. Consider the system (1). The first Lyapunov coefficient at E1 for parameter values (a, b, m, n) ∈ {(a, b, m, n)|a > 0, b > 0, m > 0, n > 1, a = a0 = b} is given by l1 =

bm3 n ln n(4b2 + 4b2 ln n + 3mn ln2 n) 4(b2 + mn ln n)(4b2 + mn ln n)2

.

(9)

Then l1 > 0, system (1) has a transversal Hopf point at E1 is unstable (weak repelling focus for the flow of system (1) restricted to the center manifold). Moreover, for each a < b, but close to b, there exists a unstable limit cycle near the stable equilibrium point E1 . Proof. From (8), the transversality condition holds at the Hopf point. Now we calculate the Lyapunov coefficient, which shows the stability of the equilibrium point and the periodic orbits which appear. Using the notation of the previous section, the multilinear symmetric functions can be written as B(x, y) = (0, m(x1 y3 + x3 y1 ), −n(1 + ln n)(x1 y2 + x2 y1 ) − n ln n(x1 y1 + x2 y2 )) ,

√

C (x, y, z ) = (0, 0, −n ln n(2 + ln n)(x1 y1 z2 + x1 y2 z1 + x1 y2 z2 + x2 y1 z1

+ x2 y1 z2 + x2 y2 z1 ) − n ln3/2 n(x1 y1 z1 + x2 y2 z2 )).

For the parameters (a, b, m, n) ∈ {(a, b, m, n)|a > 0, b > 0, m > 0, n > 1, a = a0 = b}, we have

λ1 = −2b,

√ λ2,3 = ± mn ln n i,

112

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

40 20

z

0 –20 –40 –60 10 5

2 1

0 y

0

–5 –10

–1

x

–2

Fig. 3. Singularly degenerate heteroclinic cycles of system (1) with (a, b, m, n) = (0.8696, 0.5, 0.756144, 1). Initial values (−0.001, 0.001, z (0)) and

(0.001, −0.001, z (0)), where z (0) ∈ {10.001, 20.001, 30.001}.

b

a

200

15 100 5 z

z

0

–5 –100 –15 –200 10

5 2.5

2

1

0 –2.5 –5

y

c

5

2 –1 –2

1

0

0 x

y

0

–5 –10

d

1000 500

x 10

–1 –2

x

–4

1 0.5

0 z

z

0

–500

–0.5

–1000

–1

–1500

–1.5 8

5 0 –5 y

–4

–2

0 x

2

4

4

4

2

0 y

0

–4 –8

–2 –4

x

Fig. 4. Singularly degenerate heteroclinic cycles of system (1) with (a, b, n) = (0.8696, 0.5, 1) and varying m: (a) m = 0.756144, (b) m = 0.0756144, (c) m = 0.00756144. (d) m = 0.000756144.

and

√ √  i n i n 1 − √ ,− √ , , 2 m 2 m 2   √ √ √ ib m b m + i m n ln n q= − √ , √ ,1 , √ √ 2b n + i n m ln n 2b n i − m ln nn   m(2b2 + 4b2 ln n + mn ln2 n) m(2b2 + 4b2 ln n + mn ln2 n) 2bm h11 = − ,− , 2 , 2n ln n(4b2 + mn ln n) 2n ln n(4b2 + mn ln n) 4b + mn ln n h20 = (h201 , h201 , h203 ), 

p=

where

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

a

b

2

113

4

2 z

z

1

0 0 –0.5 1

–2 2 1 0 x –1

c

1

0

x

y

–1

2

1

0

0 –1

d

3

y

–1

4

2

2 0

z

z

1 0 –1

–2

–2 2

–4 2 2 0 –2

2

0 y

x –2

x

0

0 –2

–2

y

e 4

z

2 0 –2 –4 1 2

0 x

0

–1

–2

–2

y

Fig. 5. Attractors created though the disappearance of singularly degenerate heteroclinic cycles. Parameters (a, b, m) = (0.8696, 0.5, 0.756144) and varying n: (a) n = 1.012, (b) n = 1.25, (c) n = 1.635, (d) n = 2.08, (e) n = 2.2.

√

h201 = − h202 =

√

bm(10b2 + 4b2 ln n − mn ln2 n + 6b mn ln n + 4b mn ln3/2 n)

√

6n ln n(b +

√

mn ln ni)(2b +

√

mn ln ni)2

m(−10b3 − 4b3 ln n + 12bmn ln n + 9bmn ln2 n − 26b2

√

√

√

,

mn ln n − 12b2

√

√

mn ln3/2 ni + 2(mn)3/2 ln5/2 ni)

6n ln n(b + mn ln ni)(2b + mn ln ni)2 √ 3/2 −6bm mn ln n − 8bm mn ln +8b2 mi + 8b2 m ln ni − 2m2 n ln2 ni = . √ √ 3 mn ln n( mn ln n − 2bi)2

√

h203

Then we can compute the following value G21 =

1

[3bm3 n3/2 ln3/2 n(4b2 + 4b2 ln n + 3mn ln2 n) + mn ln n)(4b2 + mn ln n)2 + im3/2 [2m3 n3 ln5 n + 6m2 n2 (3b2 + mn) ln4 n + b2 mn(48b2 + 17mn) ln3 n + 4(8b6 + 6b4 mn + 3b2 m2 n2 ) ln2 n + 16(b6 + 2b4 mn) ln n + 32b6 ]]. √

6 n ln n(

b2

Therefore, l1 =

1 2

Re G21 =

bm3 n ln n(4b2 + 4b2 ln n + 3mn ln2 n) 4(b2 + mn ln n)(4b2 + mn ln n)2

.

It is easy to find that the first Lyapunov coefficient is positive. Therefore, the theorem is proved.



,

114

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

a

0.4

Lyapunov exponents

0.2

0

–0.2

–0.4

L1 L2

–0.6

2

4

6

8

10

12

13.5

8

10

12

13.5

8

10

12

n

b

2.3

Lyapunov dimension

2

1.5

1

0.5

0

2

4

6 n

c

18 16 14 12

|s|

10 8 6 4 2 0

2

4

6 n

Fig. 6. Lyapunov exponents spectrum and bifurcation diagram of system (1) with (a, b, m) = (0.8696, 2.1, 0.756144) and 2 < n < 13.5 at the initial values (−0.08, −1.02, −0.49): (a) Lyapunov exponents; (b) Lyapunov dimension; (c) Bifurcation diagram.

4. Dynamical structure of the new chaotic system The basic dynamics of the new chaotic system can be summarized in the following singularly degenerate heteroclinic cycles, Lyapunov exponents spectrum, Lyapunov dimensions, bifurcation diagrams, and so on. 4.1. Singularly degenerate heteroclinic cycles, periodic solutions, and chaotic attractors For n = 1, the system (1) becomes x˙ = a(y − x) y˙ = −by + mxz z˙ = 1 − exy ,



(10)

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

a

b

15

15

10

10

5

5 z

z

115

0

0

–5

–5

–10 2

–10 2 0

0 y

–4

0 y

–1 x

–2 –2

c

0 –2 –4

–2

–1 x

40

z

20 0 –20 –40 5 2 0

0

y –5

–2

x

Fig. 7. Phase portraits of the system (1) with (a, m, n) = (0.8696, 0.756144, 10.5) with the initial values (−0.08,−1.02, −0.49): (a) n = 9; (b) u = 9.41; (c) u = 9.43.

which has the line of equilibria (0, 0, z ), z ∈ R. Next we are considering parameters a > 0, b > 0 and m > 0. The Jacobian matrix of system (1) at the equilibrium point (0, 0, z ) is given by

−a J = −mz 

0

a −b 0

0 0 . 0



(11)

The characteristic equation of J at (0, 0, z ) is given

λ3 + (a + b)λ2 + (ab − amz )λ = 0.

(12)

Therefore, the eigenvalues

λ1,2 =

1 2

−a − b ±



 (a − b)2 + 4amz ,

λ3 = 0,

(13)

with corresponding eigenvectors given by

v1,2 =

  −2a(a + mz ± (a − b)2 + 4amz ), −a(b − 3mz )     + (b + mz ) b + (a − b)2 + 4amz , 2(−2a + b + mz ) ,

v3 = (0, 0, 1). (a−b)2

If z < − 4am , the eigenvalues λ1,2 are complex with the negative real part. Considering the corresponding eigenvectors, this implies that the solutions locally spiraling toward the equilibrium point Q = (0, 0, z ) on a surface tangent to the plane spanned by the eigenvectors v1,2 , hence in a direction normal to the z-axis. If −

(a−b)2

< z < mb , then the eigenvalues λ1,2 real and negative. Hence, trajectories move toward the z-axis without spiraling. For z > mb , the eigenvalues λ1,2 are real with opposite signs. Then taking into account the eigenvectors v1,2 , the system has a normally hyperbolic saddle at the point P = (0, 0, z ). In the specific case in which z = mb the equilibrium point (0, 0, z ) is more degenerated, having two vanishing 4am

eigenvalues. By above analysis and detailed numerical study, when a > 0, b > 0 and m > 0, the system presents an infinite set of singularly degenerate heteroclinic cycles, which consist of invariant sets formed by a line of equilibria together with heteroclinic orbits connecting two of the equilibria. Each one of these cycles is formed by one of the one-dimensional (a−b)2

unstable manifolds of the saddle P, z > mb , which connects P with the normally hyperbolic focus Q , z < − 4am , as t → ∞. As the system presents an infinite number of normally hyperbolic saddles P and foci Q , there exists an infinite

116

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

a 1.5 1

x max

0.5 0 –0.5 –1 –1.5

–0.2

0.9

2

3.1

4.2

b

b

2 1

Lyapunov exponents

0 –1 –2 –3 –4 –5 –6 –0.2

0.8696

2 b

3.1

4.2

Fig. 8. Lyapunov exponents spectrum and bifurcation diagram of system (1) with (a, m, n) = (0.8696, 0.756144, 47) and 0.2 < b < 4.2 at initial values (−0.08, 0.02, −0.49): (a) Bifurcation diagram; (b) Lyapunov exponents.

set of singularly degenerate heteroclinic cycles. In Fig. 3, some of them are shown: for each initial condition considered sufficiently close to the saddle P at the z-axis, a singularly degenerate heteroclinic cycle is created. We also observe that the saddles P and the stable foci Q extend to infinity on the negative and positive parts of the z-axis (see Fig. 4). In the following, fixing (a, b, m) = (0.8696, 0.5, 0.756144), we consider behaviors of system (1) for the parameter n near 1. In Fig. 5(a–e), and the phase portraits are shown for n = 1.02, 1.25, 1.635, 2.08 and 2.2. These figures suggest the existence of periodic solutions and chaotic attractors when these cycles disappear. Behavior of the system evolves from periodic solution to chaotic attractor with n increasing. For example n = 2.08, we have the following Lyapunov exponents: λL1 = 0.0104, λL2 = −0.0005, and λL3 = −1.3795. 4.2. n increasing when other parameters are fixed When parameters (a, b, m) = (0.8696, 2.1, 0.756144) and varying n ∈ (2, 13.5), according to the characteristic polynomial (3) and condition (4), E1,2 both are stable. However, this system has been found to be chaotic over a wide range of parameters and has many interesting complex dynamical behaviors by varying parameter n. By the detailed numerical simulation as well as theoretical analysis, the Lyapunov exponents spectrum with initial values (−0.08, −1.02, −0.49) are shown in Fig. 6(a). The corresponding Lyapunov dimensions and bifurcation diagram are also displayed in Fig. 6(b), (c), respectively. Obviously, the maximum Lyapunov exponent is negative when n ∈ (2, 9.419), implying that the new system (1) is attracted into a sink, as shown in Fig. 7(a). When n passes through 9.419, topology structure of the system changes dramatically and the maximum Lyapunov exponent rapidly becomes positive. When n = 9.41 and n = 9.43, the phase portraits are shown in Fig. 7(b) and (c), respectively. As n increases further in the region n > 9.419, the system locates in the chaotic range.

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

a

b

100

100 50

50

0

z

0

z

117

–50

–50

–100

–100 –150 20

–150 20 2 0 –20

x

–2

0

y –20

d

150

x

–2

100

75

50

0

0

z

z

c

2 0

0

y

E1

–75

–50

–150 10

–100 10 2

2 0 y –10

e

0

0

y

x

–2

–10

f

8

8

4

4

0

0

–4

–4

–8

–8

0

50

100

150

200

–12

250

0

50

100

t

g

x

–2

12

y

y

12

–12

0

150

200

250

150

200

250

t

h

10

8 6

5

4 2 y

0

0 –2

–5

–4 –6

–10

0

50

100

150 t

200

250

300

–8 0

50

100 t

Fig. 9. Phase portraits and time series of system (1) with (a, m, n) = (0.8696, 0.756144, 47) at the initial values (−0.08, 0.02, −0.49): (a) b = 0.5; (b) b = 0.8696; (c) b = 3.2; (d) b = 3.8; (e) The time series of y(t) when b = 0.5; (f) The time series of y(t) when b = 0.8696; (g) The time series of y(t) when b = 3.2; (h) The time series of y(t) when b = 3.8.

4.3. b increasing when other parameters are fixed Now we fix (a, m, n) = (0.8696, 0.756144, 47) and vary b ∈ [−0.2,4.2]. According to the characteristic polynomial (3) and condition (4), E1,2 both are unstable for b ∈ [−0.2, 0.8696), non-hyperbolic for b = 0.8696 and stable for

118

Z. Wei, Q. Yang / Nonlinear Analysis: Real World Applications 12 (2011) 106–118

b ∈ (0.8696, 4.2]. Fig. 8(a) shows the bifurcation diagram of system (1) with respect to the parameter region. Fig. 8(b) shows the corresponding the Lyapunov exponents spectrum with initial values (−0.08, 0.02, −0.49). It is presented that two different type of chaotic attractors coexisting for certain parameter conditions. When b ∈ [−0.2, 0.8696), the max Lyapunov exponents L1 > 0 and the novel system has chaotic state with two saddle–foci. For example, with b = 0.5, the phase portrait and time series of y(t ) are shown in Fig. 9(a), (e). For b = 0.8696, system (1) is also chaotic with a positive Lyapunov exponent (see Fig. 8(b)) and corresponding phase portrait and time series of y(t) are shown in Fig. 9(b), (f). Note that for any b > 0.8696, E1,2 both are asymptotically stable, but the chaotic attractor can also be obtained in the range (0.8696, 3.75]. For example, for b = 3.2, the phase portrait and time series of y(t) are shown in Fig. 9(c), (g), respectively. When b > 3.75, as further increases, there are three negative Lyapunov exponents and the new system (1) is attracted into a sink, as shown in Fig. 9(d) and (h). 5. Conclusion In this paper, we have constructed a new chaotic system with a nonlinear term in the form of exponential function. This new chaotic attractors coexisting with stable or non-hyperbolic equilibria are different from the other attractors emerging in the Lorenz and Lorenz-like systems. Basic properties of the system have been analyzed by means of the stability and bifurcations, Lyapunov exponents spectrum, bifurcation diagram and associated Poincaré map. The existence of singularly degenerate heteroclinic cycles for a suitable choice of the parameters is investigated. Moreover, periodic solutions and chaotic attractors can be found when these cycles disappear. There are still abundant and complex dynamical behaviors and the topological structure of the new system should be completely and thoroughly investigated and exploited. Acknowledgements The authors acknowledge the referees and the editor for carefully reading this paper and suggesting many helpful comments. This work was supported by the National Natural Science Foundation of China (No. 10871074). 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]

E.N. Lorenz, Deterministic non-periodic flow, J. Atmospheric Sci. 20 (1963) 130–141. O.E. Rössler, An equation for continuous chaos, Phys. Lett. A 57 (1976) 397–398. J.C. Sprott, Some simple chaotic flows, Phys. Rev. E 50 (1994) 647–650. J.C. Sprott, A new class of chaotic circuit, Phys. Lett. A 266 (2000) 19–23. J.C. Sprott, Simplest dissipative chaotic flow, Phys. Lett. A 228 (1997) 271–274. G.R. Chen, T. Ueta, Yet another chaotic attractor, Internat. J. Bifur. Chaos 9 (1999) 1465–1466. J.H. Lü, G.R. Chen, A new chaotic attractor conined, Internat. J. Bifur. Chaos 12 (2002) 659–661. Q.G. Yang, G.R. Chen, K.F. Huang, Chaotic attractors of the conjugate Lorenz-type system, Internat. J. Bifur. Chaos 17 (2007) 3929–3949. G. van der Schrier, L.R.M. Maas, The diffusionless Lorenz equations: Sil’nikov bifurcations and reduction to an explicit map, Physica D 141 (2000) 19–36. R. Shaw, Strange attractor, chaotic behaviour and information flow, Z. Naturforsch. 36A (1981) 80–112. Q.G. Yang, G.R. Chen, A chaotic system with one saddle and two stable node-foci, Internat. J. Bifur. Chaos 18 (2008) 1393–1414. F.S. Dias, L.F. Mello, J.G. Zhang, Nonlinear analysis in a Lorenz-like system, Nonlinear Anal. RWA (2009) doi:10.1016/j.nonrwa.2009.12.010. C. Sparrow, The Lorenz Equations: Bifurcation, Chaos, and Strange Attractor, Springer-Verlag, New York, 1982. T.S. Zhou, G.R. Chen, Y. Tang, Complex dynamical behaviors of the chaotic Chen’s system, Internat. J. Bifur. Chaos 13 (2003) 2561–2574. Q.G. Yang, G.R. Chen, T.S. Zhou, A unified Lorenz-type system and its canonical form, Internat. J. Bifur. Chaos 16 (2006) 2855–2871. H. Kokubu, R. Roussarie, Existence of a singularly degenerate heterclinic cycle in the Lorenz system and its dynamical consequences: part 1*, J. Dyn. Differ. Equ. 16 (2004) 513–557. M. Messias, Dynamics at infinity and the existence of singularly degenerate heteroclinic cycles in the Lorenz system, J. Phys. A: Math. Theor. 42 (2009) 115101. L.F. Mello, S.F. Coelho, Degenerate Hopf bifurcations in the Lü system, Phys. Lett. A 373 (2009) 1116–1120. L.F. Mello, M. Messias, D.C. Braga, Bifurcation analysis of a new Lorenz-like chaotic system, Chaos Solitons Fractals 37 (2008) 1244–1255. J. Li, J. Zhang, New treatment on bifurcation of periodic solutions and homoclinic orbits at high r in the Lorenz equations, SIAM J. Appl. Math. 53 (1993) 1059–1071. Y. Yu, S. Zhang, Hopf bifurcation analysis in the Lü system, Chaos Solitons Fractals 21 (2004) 1215–1220. D. Huang, Periodic orbits and bomoclinic orbits of the diffusionless Lorenz equations, Phys. Lett. A 309 (2003) 248–253. Z.C. Wei, Q.G. Yang, Controlling the diffusionless Lorenz equations with periodic parametric perturbation, Comput. Math. Appl. 58 (2009) 1979–1987. L.P. Sil’nikov, A case of the existence of a countable number of periodic motions, Sov. Math. Docklady 6 (1965) 163–166. L.P. Sil’nikov, A contribution of the problem of the structure of an extended neighborhood of rough equilibrium state of saddle–focus type, Math. USSR-Shornik 10 (1970) 91–102. C.P. Silva, Sil’nikov theorem-a tutorial, IEEE Trans. Circuits Syst. I 40 (1993) 657–682. Y. Hardy, W.H. Steeb, The Rikitake two-disk dynamo system and domains with periodic orbits, Int. J. Theor. Phys. 38 (1999) 2413–2417. J. Llibre, M. Messias, Global dynamics of Rikitake system, Physica D 238 (2009) 241–252. J. Llibre, X. Zhang, Invariant algebraic surfaces of the Rikitake system, J. Phys. A: Math. Gen. 33 (2000) 7613–7635.