Analysis of a New Quadratic 3D Chaotic Attractor

Abstract

A new three-dimensional chaotic system is introduced. Basic properties of this system show that its corresponding attractor is topologically different from some well-known systems. Next, detailed information on dynamic of this system is obtained numerically by means of Lyapunov exponents spectrum, bifurcation diagrams, and 0-1 chaos indicator test. We finally prove existence of this chaotic attractor theoretically using Shil’nikov theorem and undetermined coefficient method.

1. Introduction

The current surge towards the study on dynamical systems and 3D chaotic attractors started by the remarkable discovery of Lorenz in 1963 [1]. While he was studying atmospheric convections using Saltzman equations, he finally came up with a new system of differential equations (1), which is known after his name. His further investigations on this new system showed that for specific values of parameters it has a new type of attractor, namely, chaotic attractor as follows: ̇,̇̇𝑥=𝑎𝑦−𝑥𝑦=𝑐𝑥−𝑦−𝑥𝑧,𝑧=𝑥𝑦−𝑏𝑧. (1)

During the past decades, enormous amount of researches have been done on this system which have revealed characteristics and features of it [2].

Besides the researches on Lorenz, many other system of equations have been introduced and analysed since then, like Chen, Rössler, Chua, Lü, Qi, and many others. For instance, a new chaotic system was reported with no equilibria which is based on Sprott D system [3]. In another research, a chaotic system with only one stable equilibrium was introduced which typically is not anticipated to show such a behaviour [4]. And recently a technique to construct a chaotic system with an arbitrary number of equilibria has been also introduced [5]. Since the idea of the proposed system has come from Chen and Qi systems, we describe these two systems briefly.

The Chen system (2) is constructed using a state feedback in the second equation of Lorenz, and he showed that this system is chaotic at 𝑎=35, 𝑏=3, and 𝑐=28 as follows: ̇,̇()̇𝑥=𝑎𝑦−𝑥𝑦=𝑐−𝑎𝑥+𝑐𝑦−𝑥𝑧,𝑧=𝑥𝑦−𝑏𝑧. (2)

Later on, by eliminating the first term in the second equation, Chen and Lü proposed the Lü system, which is infact the transition between Lü and Chen [6].

In 2005, Qi and his colleagues added a cross-product nonlinear term to the first equation of Lorenz and introduced a new system (3) that was topologically different from Lorenz, Chen, Rössler, and even Lorenz system family [7] as follows: ̇̇̇𝑥=𝑎𝑦−𝑥+𝑦𝑧,𝑦=𝑐𝑥−𝑦−𝑥𝑧,𝑧=𝑥𝑦−𝑏𝑧. (3) Comparing (2) and (3) with Lorenz shows that how small modifications on Lorenz could lead to a system with new characteristic features and in some cases even topologically different attractor. The system we are going to introduce and study here is ̇̇()̇𝑥=𝑎𝑦−𝑥+𝑦𝑧,𝑦=𝑐−𝑎𝑥+𝑐𝑦−𝑥𝑧,𝑧=𝑥𝑦−𝑏𝑧. (4)

It has a chaotic attractor at 𝑎=15, 𝑏=8/3, and 𝑐=10 which is topologically different from Lorenz, Chen, and even Rössler to be shown shortly.

This paper is divided into five parts. Section 2 deals with studying the basic properties of system (4). After that in Section 3, more detailed investigation will be done on dynamic of this system using Lyapunov exponents spectrum, bifurcation diagrams, and 0-1 test which reveals the behaviour of this system for various sets of parameters, and how it evolves from one state to the other. Then in Section 4 we will theoretically prove the existence of chaos in this system using Shil'nikov theorem and undetermined coefficient method [8]. We conclude this paper in Section 5.

2. Basic Properties

In this section, we start with equilibrium points of the system and check their stability at initial values of parameters 𝑎, 𝑏, and 𝑐. Putting equations of the system equal to zero gives the equilibrium points 𝑣₁=(),𝑣0,0,0_2,3=±𝑎√2√𝑐𝑑𝑏𝑎((√2𝑐𝑎−𝑑𝑎/2𝑐)),±2√𝑐𝑑𝑏𝑎,𝑎2𝑐²𝑑((,𝑣2𝑐𝑎−𝑑𝑎/2𝑐))_4,5=±𝑎√−2𝑐𝑒𝑏𝑎((√2𝑐𝑎+𝑒𝑎/2𝑐)),±−2𝑐𝑒𝑏𝑎𝑎2𝑐,−²𝑒((,2𝑐𝑎+𝑒𝑎/2𝑐)) (5) where 𝑑=3𝑐+√𝑐²√+4𝑐𝑎,𝑒=−3𝑐+𝑐²+4𝑐𝑎. (6) Therefore, at the initial values of parameters 𝑎, 𝑏, and 𝑐, this system has five equilibrium points in contrast to Lorenz, Chen, and Rössler which have three, three, and two equilibriums, respectively. It implies topological distinction between system (4) and these systems.

In order to check their stability we derive Jacobian matrix of the system 𝐷𝑓𝑥,𝑦,𝑧=−𝑎𝑎+𝑧𝑦𝑐−𝑎−𝑧𝑐−𝑥𝑦𝑥−𝑏. (7) Substituting the origin in this matrix we get the characteristic equation 𝜆²+(−𝑐+𝑎)𝜆+𝑎(−2𝑐+𝑎)(𝜆+𝑏)=0 (8) which has three eigenvalues 𝜆₁=6.5139, 𝜆₂=−2.6667, and 𝜆₃=−11.5139. Since all the eigenvalues are real, Hartman-Grobman theorem implies that origin is a saddle point which is not Lyapunov stable according to the Lyapunov theorem of stability.

The eigenvalues of the Jacobian at 𝑣₂ and 𝑣₃ are 𝜆₁=−19.5735 and 𝜆_2,3=5.9534±11.4898𝑖; notice that these two points are symmetrically located around 𝑧 axis. None of these numbers have real part zero, and because 𝜆_2,3 are complex, 𝑣₂ is an unstable saddle-focus point. Exactly the same results holds for 𝑣₃.

At 𝑣₄ and 𝑣₅ the eigenvalues are 𝜆₁=−10.3231 and 𝜆_2,3=1.3282±6.5092𝑖. Similar to 𝑣_2,3, it shows that this system has a saddle focus at 𝑣₄ and 𝑣₅ that is not Lyapunov stable.

A large class of chaotic systems are dissipative dynamical systems which satisfy the condition 𝜕𝑥̇∇⋅𝑉=+𝜕𝑦̇𝜕𝑥+𝜕𝑧̇𝜕𝑦𝜕𝑧<0. (9) In the case of system (4) we have ∇⋅𝑉=𝑐−(𝑎+𝑏)=10−(15+8/3)<0; therefore dissipativity condition holds on this system. Moreover 𝑑𝑉𝑑𝑡=𝑒^{𝑐−(𝑎+𝑏)}=0.0004682. (10) It implies that the volume of the attractor decreases by a factor of 0.0004682 and at each particular time 𝑡 the volume is 𝑉(𝑡)=𝑉₀𝑒^{(𝑐−(𝑎+𝑏))𝑡}. (11)

The key characteristic of a chaotic system is its sensitive dependency on initial conditions; which is associated with having positive Lyapunov exponents. Using the MATLAB program at [9], which is based on the algorithm introduced in [10], the Lyapunov exponents of the system are 0.766±0.005, 0±0.0004, and −8.435±0.005. In addition, exponentiating these numbers indicates that nearby trajectories diverge from each other by a factor of 2.151.

Having the Lyapunov exponents, Kaplan-Yorke (Lyapunov) dimension [11] of the attractor is 2.0908. On the other hand, using correlation and box-counting methods, fractal dimensions of the attractor are 1.9763 and 0.9405, respectively, which indicate its fractal structure. Figure 1 shows the attractor and its projections on different planes.

(a)

(b)

(c)

(d)

Figure 1 

3D view of the attractor and its various projections.

3. Computational Analysis

In this section, we employ extensive computations to have a broader picture of the dynamic of this system in different regions of parameters' values. Our focus is on the parameter 𝑎, but we do show the results regarding 𝑏 and 𝑐 in the figures.

At different set of parameters' values, behaviour of a system could be either chaotic, periodic, or convergence to one of the equilibria. So to have a good understanding of its dynamics, we need to identify the regions whereby the system is in either of these states. Here we use Lyapunov exponents spectrum, bifurcation diagrams, and 0-1 chaos indicator test [12] to obtain this information.

Figure 2 shows the Lyapunov exponents spectrum with respect to parameter 𝑎∈[11.5,25]. As we see, at the beginning the system is in periodic state, and trajectories converge to a limit cycle in the state space. At 𝑎=13.94, chaotic state emerges, which means that the system becomes sensitive to its initial condition. It continues up to 17.68 where the maximum exponent becomes negative and trajectories fall to convergence to the equilibria.

Figure 2 

Lyapunov exponents spectrum of the system with respect to 𝑎.

This result can be justified using bifurcation diagram and 0-1 test. Figure 3 shows that following the periodic state, the system experiences a period doubling cascade leading to chaos right in the interval 𝑎∈[13.94,17.68]. The window in this region corresponds to the small region in the Lyapunov spectrum whereby the maximum exponents become zero temporarily. Figure 4 shows the result of 0-1 test which is clearly in accordance with Figures 2 and 3. Results regarding parameters 𝑏 and 𝑐 are displayed in Figures 5, 6, 7, and 8.

Figure 3 

Bifurcation diagram of the system with respect to 𝑎.

Figure 4 

Implication of 0-1 test on the system with respect to 𝑎.

Figure 5 

Lyapunov exponents spectrum of the system with respect to 𝑏.

Figure 6 

Bifurcation diagram of the system with respect to 𝑏.

Figure 7 

Lyapunov exponents spectrum of the system with respect to 𝑐.

Figure 8 

Bifurcation diagram of the system with respect to 𝑐.

4. Theoretical Proof of Chaos

4.1. Finding the Heteroclinic Orbit

Numerical simulations do not give us the rigorous proof of the existence of chaos in a system. Therefore, we have to work theoretically to find the final answer to this question. The method we go through here is Shil'nikov theorem for the existence of chaos. Having the conditions of this theorem satisfied, it says that the system has Smale horseshoes and horseshoe type of chaos.

Consider the three-dimensional autonomous system 𝑥̇=𝑓(𝑥), (12) where vector field 𝑓∶ℝ³→ℝ³ belongs to class ℂ^𝑟(𝑟≥1). Shil'nikov theorem is [13].

Theorem 1 (Heteroclinic Shil'nikov theorem). Suppose that two distinct equilibrium points, denoted by 𝑃₁ and 𝑃₂, respectively, of system 𝑥̇=𝑓(𝑥), are saddle foci whose characteristic values 𝛾_𝑘, 𝜌_𝑘+𝑖𝜔_𝑘, (𝑘=1,2) satisfy the Shil'nikov inequalities 𝜌₁𝜌₂>0 or 𝛾₁𝛾₂>0. Suppose also that there exists a heteroclinic orbit jointing 𝑃₁ and 𝑃₂. Then the system 𝑥̇=𝑓(𝑥) has both Smale horseshoes and the horseshoe type of chaos.

Since the origin is not a saddle focus, Shil'nikov theorem cannot be applied on it. However, 𝑣_2,3 and 𝑣_4,5 satisfy the theorem's condition. Here we put 𝑣_4,5 in our consideration and show the presence of a heteroclinic orbit between them which proves the existence of chaos.

The method we employ to find the heteroclinic orbit is undetermined coefficients method [8] with a new approach [14]. Without loss of generality, assume the direction from 𝑣₄ to 𝑣₅ corresponding to 𝑡→∞ and from 𝑣₅ to 𝑣₄ corresponding to 𝑡→−∞.

In the case of 𝑡>0 and according to the undetermined coefficients method, assume that 𝑥(𝑡)=𝑥₅+^∞_𝑘=1𝑎_𝑘𝑒^𝑘𝛼𝑡,𝑦(𝑡)=𝑦₅+^∞_𝑘=1𝑏_𝑘𝑒^𝑘𝛼𝑡,𝑧(𝑡)=𝑧₅+^∞_𝑘=1𝑐_𝑘𝑒^𝑘𝛼𝑡,𝑡>0. (13) Substituting these equations into the system (4) we have ^∞_𝑘=1𝑎_𝑘𝑘𝛼𝑒^𝑘𝛼𝑡=𝑎𝑦₅+^∞_𝑘=1𝑏_𝑘𝑒^𝑘𝛼𝑡−𝑥₅−^∞_𝑘=1𝑎_𝑘𝑒^𝑘𝛼𝑡+𝑦₅𝑧₅+𝑦₅^∞_𝑘=1𝑐_𝑘𝑒^𝑘𝛼𝑡+𝑧₅^∞_𝑘=1𝑏_𝑘𝑒^𝑘𝛼𝑡+^∞_{𝑖+𝑗=𝑘≥2}𝑏_𝑖𝑐_𝑗𝑒^𝑘𝛼𝑡,_∞𝑘=1𝑏_𝑘𝑘𝛼𝑒^𝑘𝛼𝑡=(𝑐−𝑎)𝑥₅+(𝑐−𝑎)^∞_𝑘=1𝑎_𝑘𝑒^𝑘𝛼𝑡+𝑐𝑦₅+𝑐^∞_𝑘=1𝑏_𝑘𝑒^𝑘𝛼𝑡−𝑥₅𝑧₅−𝑥₅^∞_𝑘=1𝑐_𝑘𝑒^𝑘𝛼𝑡−𝑧₅^∞_𝑘=1𝑎_𝑘𝑒^𝑘𝛼𝑡−^∞_{𝑖+𝑗=𝑘≥2}𝑎_𝑖𝑐_𝑗𝑒^𝑘𝛼𝑡,_∞𝑘=1𝑐_𝑘𝑘𝛼𝑒^𝑘𝛼𝑡=𝑥₅𝑦₅+𝑥₅^∞_𝑘=1𝑏_𝑘𝑒^𝑘𝛼𝑡+𝑦₅^∞_𝑘=1𝑎_𝑘𝑒^𝑘𝛼𝑡+^∞_{𝑖+𝑗=𝑘≥2}𝑎_𝑖𝑏_𝑗𝑒^𝑘𝛼𝑡−𝑏𝑧₅−𝑏^∞_𝑘=1𝑐_𝑘𝑒^𝑘𝛼𝑡. (14)

Comparing the coefficients of 𝑒^𝑘𝛼𝑡 of the same power we have 𝛼𝐼−𝐷𝑓𝑣₅⋅𝑎₁𝑏₁𝑐₁𝑣=0,for𝑘=1,𝑘𝛼𝐼−𝐷𝑓₅⋅𝑎_𝑘𝑏_𝑘𝑐_𝑘=_𝑖+∑_𝑗=𝑘𝑏_𝑖𝑐_𝑗−∑_{𝑖+𝑗=𝑘}𝑎_𝑖𝑐_𝑗𝑖+∑_𝑗=𝑘𝑎_𝑖𝑏_𝑗 (15)(16) which holds for 𝑘≥2. In (15), (𝑎₁,𝑏₁,𝑐₁)≠0 because otherwise, by induction we get 𝑎_𝑘,𝑏_𝑘,𝑐_𝑘=0 (17) for all 𝑘≥2. In addition, since jacobian at 𝑣₅ has a negative eigenvalue, there is a unique 𝛼 such that det𝛼𝐼−𝐷𝑓𝑣₅=0. (18) Thus, we can easily determine 𝑎₁, 𝑏₁, and 𝑐₁ using a free variable 𝜃. Also note that in (16) det𝑘𝛼𝐼−𝐷𝑓𝑣₅≠0. (19)

Now we consider 𝑡<0, corresponding to moving from 𝑣₅ to 𝑣₄. Using the variable transformation 𝜏=−𝑡 with 𝑡>0, system (4) transforms to 𝑑𝑥𝑑𝜏=𝑎𝑥−𝑦−𝑦𝑧,𝑑𝑦=()𝑑𝜏𝑎−𝑐𝑥−𝑐𝑦+𝑥𝑧,𝑑𝑧𝑑𝜏=−𝑥𝑦+𝑏𝑧. (20) Similarly assume that 𝑥(𝑡)=𝑥₄+^∞_𝑘=1𝑎_𝑘𝑒^{−𝑘𝛽𝜏},𝑦(𝑡)=𝑦₄+^∞_𝑘=1𝑏̂_𝑘𝑒^{−𝑘𝛽𝜏},𝑧(𝑡)=𝑧₄+^∞_𝑘=1𝑐_𝑘̂𝑒^{−𝑘𝛽𝜏},𝜏<0. (21) Substituting these equations into the system (20) we obtain ^∞_𝑘=1−𝑎_𝑘𝑘𝛽𝑒^{−𝑘𝛽𝜏}=𝑎𝑥₄+𝑎𝑎_𝑘𝑒^{−𝑘𝛽𝜏}−𝑎𝑦₄−𝑎^∞_𝑘=1𝑏̂_𝑘𝑒^{−𝑘𝛽𝜏}−𝑦₄𝑧₄−𝑦₄^∞_𝑘=1𝑐_𝑘̂𝑒^{−𝑘𝛽𝜏}−𝑧₄^∞_𝑘=1𝑏̂_𝑘𝑒^{−𝑘𝛽𝜏}−^∞_{𝑖+𝑗=𝑘≥2}𝑐_𝑖̂𝑏̂_𝑗𝑒^{−𝑘𝛽𝜏},_∞𝑘=1−𝑏̂_𝑘𝑘𝛽𝑒^{−𝑘𝛽𝜏}=(𝑎−𝑐)𝑥₄+^∞_𝑘=1𝑎_𝑘𝑒^{−𝑘𝛽𝜏}−𝑐𝑦₄−^∞_𝑘=1𝑏̂_𝑘𝑒^{−𝑘𝛽𝜏}+𝑥₄𝑧₄+^∞_𝑘=1𝑐_𝑘̂𝑒^{−𝑘𝛽𝜏}+𝑧₄^∞_𝑘=1𝑎_𝑘𝑒^{−𝑘𝛽𝜏}+^∞_{𝑖+𝑗=𝑘≥2}𝑒^{−𝑘𝛽𝜏},_∞𝑘=1−𝑐_𝑘̂𝑘𝛽𝑒^{−𝑘𝛽𝜏}=−𝑥₄𝑦₄−𝑥₄^∞_𝑘=1𝑏̂_𝑘𝑒^{−𝑘𝛽𝜏}−𝑦₄^∞_𝑘=1𝑎_𝑘𝑒^{−𝑘𝛽𝜏}−^∞_{𝑖+𝑗=𝑘≥2}𝑎_𝑖𝑏̂_𝑗𝑒^{−𝑘𝛽𝜏}+𝑏𝑧₄+𝑏^∞_𝑘=1𝑐_𝑘̂𝑒^{−𝑘𝛽𝜏}. (22)

Now we equate the coefficients of 𝑒^{−𝑘𝛽𝜏} which have same power terms. It gives us the following: 𝛽𝐼−𝐷𝑓𝑣₄⋅𝑎₁𝑏̂₁𝑐₁̂𝑣=0,for𝑘=1,𝑘𝛽𝐼−𝐷𝑓₄⋅𝑎_𝑘𝑏̂_𝑘𝑐_𝑘̂=−∑_{𝑖+𝑗=𝑘}𝑏̂_𝑖𝑐_𝑗̂_𝑖+∑_𝑗=𝑘𝑎_𝑖𝑐_𝑗̂−∑_{𝑖+𝑗=𝑘}𝑎_𝑖𝑏̂_𝑗,for𝑘≥2. (23) Similar to the previous case, 𝑎₁, 𝑏̂₁, and 𝑐₁̂ can be determined by a free variable and det(𝑘𝛽𝐼−𝐷𝑓(𝑣₄))≠0.

So the final expression of the heteroclinic orbit connecting 𝑣₄ and 𝑣₅ is 𝑥(𝑡)=𝑥₅+∑^∞_𝑘=1𝑎_𝑘(𝑎,𝑏,𝑐,𝛼)𝑒^𝑘𝛼𝑡𝑥,𝑡>0,₄+∑^∞_𝑘=1𝑎_𝑘𝑎,𝑏,𝑐,𝛽𝑒^{−𝑘𝛽𝑡}𝑦(𝑡)=𝑦,𝑡<0,₅+∑^∞_𝑘=1𝑏_𝑘(𝑎,𝑏,𝑐,𝛼)𝑒^𝑘𝛼𝑡𝑦,𝑡>0,₄+∑^∞_𝑘=1𝑏̂_𝑘𝑎,𝑏,𝑐,𝛽𝑒^{−𝑘𝛽𝑡}𝑧(𝑡)=𝑧,𝑡<0,₅+∑^∞_𝑘=1𝑐_𝑘(𝑎,𝑏,𝑐,𝛼)𝑒^𝑘𝛼𝑡𝑧,𝑡>0,₄+∑^∞_𝑘=1𝑐_𝑘̂𝑎,𝑏,𝑐,𝛽𝑒^{−𝑘𝛽𝑡},𝑡<0. (24)

4.2. Convergence of the Series Solution

Since the solution is in series form, we have to prove that it is convergence before making any conclusion. Having (15) and (16) and the fact that (𝑘𝛼𝐼−𝐷𝑓(𝑣₅)) is nonsingular, for any 𝑘≥1 we have 𝑎_𝑘𝑏_𝑘𝑐_𝑘=𝑚_𝑎𝜃^𝑘𝑚_𝑏𝜃^𝑘𝑚_𝑐𝜃^𝑘, (25) where 𝑚_𝑎, 𝑚_𝑏, and 𝑚_𝑐 are constants and 𝜃 is the free variable we use to determine 𝑎₁, 𝑏₁, and 𝑐₁. In addition, after doing some algebra on (16) using the default values of 𝑎, 𝑏, and 𝑐 we have 𝑎_𝑘<𝑀𝜃^𝑘𝑘,𝑘≥2. (26)

As we mentioned before, 𝜃 is a free variable, so we can assume it to be 0<𝜃<1. Thus 𝑥(𝑡)=𝑥₅+^∞_𝑘=1𝑎_𝑘𝑒^𝑘𝛼𝑡<𝑥₅+𝑎₁𝑒^𝛼𝑡+^∞_𝑘=2𝑀𝜃^𝑘𝑘𝑒^𝑘𝛼𝑡<𝑥₅+𝑎₁𝑒^𝛼𝑡+𝑀𝜃^∞_𝑘=2𝑒^𝑘𝛼𝑡𝑘=𝑥₅+𝑎₁𝑒^𝛼𝑡+𝑀𝜃−ln1−𝑒^𝛼𝑡−𝑒^𝛼𝑡<𝑥₅+𝑎₁+𝑀𝜃. (27) Exactly the same result holds for the second part of 𝑥(𝑡) corresponding to 𝑡<0. Therefore the convergence of the series solution of 𝑥(𝑡) is proved. Similar procedure can be done to prove the convergence of 𝑦(𝑡) and 𝑧(𝑡).

One remark is necessary to give here. If we consider 𝑣₂ and 𝑣₃, we see that they satisfy in Shil'nikov theorem condition and we can follow just the same algebraic procedure to find the heteroclinic orbit between them. However, the location of these two equilibrium points are very far from the chaotic attractor observed. Therefore more investigation is required to find out the role of this heteroclinic orbit on the dynamic of the chaotic system (4).

5. Conclusion

In this paper, a new 3D quadratic chaotic attractor was introduced and analysed. We studied basic properties of this system by means of Jacobian matrices, Lyapunov exponents, and various definitions of fractal dimension which showed topological distinction between this attractor and some other well-known systems. In addition, more detailed information on its dynamic was obtained using Lyapunov exponents spectrum, bifurcation diagrams, and 0-1 chaos indicator test that disclosed just a portion of the underlying dynamic of this system. And finally, we proved the existence of chaos by showing that this system has Smale horseshoes and horseshoe type of chaos.