Bifurcation and chaos in discrete FitzHugh–Nagumo system

Bifurcation and chaos in discrete FitzHugh–Nagumo system

Chaos, Solitons and Fractals 21 (2004) 701–720 www.elsevier.com/locate/chaos Bifurcation and chaos in discrete FitzHugh–Nagumo system q Zhujun Jing b...

763KB Sizes 0 Downloads 65 Views

Chaos, Solitons and Fractals 21 (2004) 701–720 www.elsevier.com/locate/chaos

Bifurcation and chaos in discrete FitzHugh–Nagumo system q Zhujun Jing b

a,b,*

, Yu Chang c, Boling Guo

c

a Department of Mathematics, Hunan Normal University, Hunan Changsha 410081, PR China Academy of Mathematics and Systems Sciences, Chinese Academy of Sciences, Beijing 100080, PR China c Institute of Applied Physics and Computational Mathematics, Beijing 100088, PR China

Accepted 15 December 2003

Abstract The discrete FitzHugh–Nagumo system obtained by Euler method is investigated. The conditions of existence for fold bifurcation, flip bifurcation and Hopf bifurcation are derived by using center manifold theorem and bifurcation theory, chaotic behavior in the sense of Marotto’s definition of chaos is proved. And numerical simulation results not only show the consistence with the theoretical analysis but also display the new and interesting dynamical behaviors, including attracting invariant circle, period-3, period-6, period-7, period-9, period-15, period-20, period-21, and periodn orbits, an inverse cascade of period-doubling bifurcation in period-3, cascade of period-doubling bifurcation in periods-9, 15, 20 and 21, interior and exterior crisis phenomena, intermittency mechanic, transient chaos in periodwindow, attracting and non-attracting chaotic attractors. The computations of Lyapunov exponents confirm the chaotic behaviors. Ó 2004 Elsevier Ltd. All rights reserved.

1. Introduction Hodgkin and Huxley [14] constructed a HH nerve model of four non-linear differential equations to describe the electrical activity that underlies the conduction of a nerve impulse in the squid giant axon. The HH model was investigated by using different methods [1,10,20–23]. Because of the complexity FitzHugh–Nagumo (FHN) [5,21] provided a simplified version to the HH model which capture key features of the full system. And FHN model can be taken as a representative of a wide class of non-linear excitable-oscillatory system and have wide applications in the modelling of biological processes which is clearly a very important field. The modelling of dynamics of the heart is a classical example, which was already studied in [27]. The modelling of cardiac excitable medium is given by using FHN equations which is studied in [30]. The other biological processes that show in their behaviors some relationship with BVP equations are the Bonhoeffer–Van der Pol (BVP) system and CAMP signaling system and the model of the cell cycle [2]. The dynamical behaviors under different parameter conditions for FHN model have been studied in the mathematics literature. The global asymptotic stability of an equilibrium point is given in [3,7]. The conditions of the non-existence of homoclinic orbits and the existence of unique non-trivial closed orbit(the FHN system is transformed to the Lienard system) are determined in [11–13]. The bifurcation manifolds concerning the connections between saddle which are including the homoclinic bifurcation, Hopf bifurcation, fold bifurcation and Bogdanov–Takens bifurcation, are proved by using bifurcation theory and qualitative analysis in [24]. The stability analysis for the travelling wave solution and global attractor of the FHN of reaction–diffusion equations are obtained in [6,7,16,26,29]. Moreover, results in [2,15]

q

This work was supported by Chinese Academy Sciences (KZCX2-SW-118). Corresponding author. Tel.: +86-10-1-6265-1304. E-mail address: [email protected] (Z. Jing).

*

0960-0779/$ - see front matter Ó 2004 Elsevier Ltd. All rights reserved. doi:10.1016/j.chaos.2003.12.043

702

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

indicate that there are complex dynamical behaviors including (quasi-)periodic orbits, ten-periodic orbits, a cascade of period-doubling bifurcation, and chaotic behavior in Marotto’s chaos and intermittence’s chaos. In our paper, we apply the forward Euler Scheme to discrete the FHN model where the model is considered as a general model for an excitable system, such as a single cardiac cell [30], and investigate this version as discrete dynamical system in R2 by using bifurcation theory of continuous and discrete systems in [8,9,18,19,25,28]. We first rigorously prove that this discrete model possesses the fold bifurcation, flip bifurcation and Hopf bifurcation and chaos in the sense of Marotto’s definition [17] and give the bifurcation diagrams. Specifically, we analysis the case of existing three fixed points for the model and give the results of the numerical simulations that not only perfectly show the consistence with the theoretical analysis but also find the new and interesting dynamical behaviors including attracting invariant circle, period-3, period-6, period-7, period-9, period-15, period-20, period-21, and period-n orbits, an inverse cascade of period-doubling bifurcation in period-3, cascade of period-doubling bifurcation in periods-9, 15, 20 and 21, interior and exterior crisis chaotic phenomena, intermittency mechanic, attracting and non-attracting chaotic attractors. Meanwhile the computations of Lyapunov exponents confirm the chaotic behaviors. This paper is organized as follows. In Section 2, we analysis existence and stability of the fixed points and give sufficient conditions of existence for fold bifurcation, flip bifurcation and Hopf bifurcation. In Section 3, we first rigorously prove the existence of chaos in the sense of Marotto’s definition. The results of numerical simulations and the computations of Lyapunov exponents are presented in Section 4 to verify the theoretical analysis and display the complex and interesting dynamics.

2. Existence and stability of fixed points and bifurcations 2.1. Existence and stability of fixed points Consider the FHN equation [5,30] which takes the form ( 3 du ¼ f ðu; vÞ ¼ uu e=3v ; dt dv ¼ gðu; vÞ ¼ eðu þ b  cvÞ; dt

ð1Þ

where 0 < e < 1, b > 0 and c > 0 are constants. The variable u is the voltage across the neural membrane, v is a quantity of refractoriness, and the system (1) can be taken as a general for an excitable system, such as a single cardiac cell [30]. Applying Euler scheme to system (1), we obtain the discrete system     x x þ dðx  x3 =3  yÞ=e ! ; ð2Þ F : y y þ deðx þ b  cyÞ where d is the integral step size and it is actually very small, so let 0 < d < 1. The fixed points of map (2) satisfy the following equations   ( x3 þ 3  1 þ 1c x þ 3bc ¼ 0; y ¼ 1c ðx þ bÞ:

ð3Þ

By a simple analysis, it is easy to obtain the following proposition: Proposition 1 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffi þ D þ 3  3b  D, (1) If 9b2 c þ 4ð1  cÞ3 > 0, then the map (2) has an unique fixed point at Zðx; yÞ, where x ¼ 3  3b 2c 2c 3 2 D ¼ 9b cþ4ð1cÞ , and y ¼ 1c ðx þ bÞ. qffiffiffiffi 4c3 2 , (2) If 9b c þ 4ð1  cÞ3 ¼q0,ffiffiffiffi then the map (2) has two fixed points at Z1 ðx1 ; y1 Þ and Z2 ðx2 ; y2 Þ, where x1 ¼ 2 3 3b 2c y1 ¼ 1c ðx1 þ bÞ, x2 ¼ 2

3

3b , 2c

y2 ¼ 1c ðx2 þ bÞ.

3

(3) If 9b c þ 4ð1  cÞ < 0, the map (2) has three fixed points. The curve C ¼ fðb; cÞj 9b2 c þ 4ð1  cÞ3 ¼ 0g in (b; c) plane is shown in Fig. 1. There are two fixed points on the curve C, and three fixed points in the region (I), and unique fixed point in region (II). We now investigate the linear stability of the fixed points for (2).

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

703

Fig. 1. Fixed points for map (2) in parameter space (b–c).

The Jacobian matrix J of the map (2) evaluated at fixed point Zðx; yÞ is given by   1 þ dð1  x20 Þ=e d=e J¼ de 1  dec

ð4Þ

and the characteristic equation of the Jacobian matrix J can be written as k2 þ pðxÞk þ qðxÞ ¼ 0;

ð5Þ

where p ¼ 2 þ dec  dð1  x2 Þ=e; q ¼ 1  dec þ dð1  x2 Þ=e  d2 cð1  x2 Þ þ d2 : A simple calculation shows the stability of fixed points as the following. Proposition 2. The fixed point Z ¼ ðx; yÞ of map (2) is stable if the x-coordinate of this fixed point satisfies one of the following conditions: 2

Þ=e (1) 1  2e þ e2 c 6 x2 6 1 þ 2e þ e2 c, 1  cð1  x2 Þ 6¼ 0, d < ecð1x ; 1cð1x2 Þ de ;  2ed  2dec g < x2 < 1  maxfeð2  ecÞ; eðdec4Þ g, 2  edc < 0; (2) 1  minfeð2þedcÞ d d de g < x2 < 1  eðdec4Þ g, 2  edc < 0; (3) 1  minfeð2 þ dcÞ;  2ed  2dec d de ; eð2 þ ecÞg < x2 < 1  maxf 2ed  2dec ; eðdec4Þ g, 2  edc > 0; (4) 1  minfeð2þedcÞ d d 2 2e (5) 2  edc ¼ 0, and 1 < x < 1  maxf d ; eð2  ecÞg; (6) 2  edc ¼ 0, and 1 þ eð2 þ ecÞ < x2 < 1 þ 2ed;

; eð2  ecÞg; (7) 1  minfe2 c; 1cg < x2 < 1  maxfeð2þdecÞ d . (8) 1 þ eð2 þ ecÞ < x2 < 1  eð2þdecÞ d 2.2. Bifurcations In analysis of fold bifurcation and flip bifurcation b is as bifurcation parameter. qffiffiffiffiffiffiffiffiffiffi 1 1  for c > 1, thus there are the eigenvalues From 9b2 c þ 4ð1  cÞ3 ¼ 0, we have b0 ¼ 2ðc1Þ c 3 k1 ¼ 1 and

k2 ¼ 1  pðb0 Þ ¼ 1  dec þ dð1  x20 Þ=e:

704

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

qffiffiffiffiffi qffiffiffiffiffi  For the fixed point Z2 ðxðb0 Þ; yðb0 ÞÞ, where xðb0 Þ ¼ x0 ¼ 3 3b2c0 , and yðb0 Þ ¼ y0 ¼ 1c 3 3b2c0 þ b0 . Next we will prove that Z2 ðxðb0 Þ; yðb0 ÞÞ is a fold bifurcation point under certain conditions by center manifold theorem and bifurcation analysis. We first require jk2 j 6¼ 1, it leads to ec 6¼ 1 and

2  dec þ

d 6¼ 0: ce

ð6Þ

~ ¼ b  b , we transform the fixed point Z2 to the origin, and consider the Let u ¼ x  xðb0 Þ, v ¼ y  yðb0 Þ and b 0 ~ as a new and dependent variable, then map (2) becomes: parameter b   1 0 0 1 u3 þ3u2 x0 þ3ux20 uþd u  v =e u 3 C ~A ! B @b ð7Þ @ A: ~ b v ~  cvÞ v þ deðu þ b Let 0

1

c

B T ¼ @0 1

1

1c2 e2 de dc2 e dc

1

0 C A ce2

and use the translation 0 1 0 1 ~x u ~ @ b A ¼ T @ l A; ~y v then the map (7) becomes 0 1 1 0 10 1 0 ~x ~x f ð~x; l; ~y Þ 1 1 0 C B B C C AB 0 0 @lA þ @ @lA ! @0 1 A; ~y ~y gð~x; l; ~y Þ 0 0 1  pðb0 Þ

ð8Þ

where f ð~x; l; ~y Þ ¼

dec ð3x0 ðc~x þ l þ ~y Þ2 þ ðc~x þ l þ ~y Þ3 Þ; 3ð1  c2 e2 Þ

gð~x; l; ~y Þ ¼  l¼

d ð3x0 ðc~x þ l þ ~y Þ2 þ ðc~x þ l þ ~y Þ3 Þ; 3eð1  c2 e2 Þ

de ~ b: 1  c2 e2

By center manifold theory, we know that the stability of ð~x; ~y Þ ¼ ð0; 0Þ near l ¼ 0 can be determined by studying a one-parameter family of equation on a center manifold, which can be represented as follows W c ð0Þ ¼ fð~x; l; ~y Þ 2 R3 j~y ¼ hð~x; lÞ; hð0; 0Þ ¼ 0; Dhð0; 0Þ ¼ 0g; for ~x and l sufficiently small. We assume a center manifold of the form hð~x; lÞ ¼ a1 l2 þ a2~xl þ a3~x2 þ Oððj~xj þ jljÞ3 Þ:

ð9Þ

The center manifold must satisfy Nðhð~x; lÞÞ ¼ hð~x þ l þ f ð~x; l; hð~x; lÞÞ; lÞ  ð1  pðb0 ÞÞhð~x; lÞ  gð~x; l; hð~x; lÞÞ ¼ 0:

ð10Þ

Since we only concern the terms with orders lower than 3 for map (8), we can obtain the map restricted to the center manifold directly from (8), which is given by dec3 x0 2 decx0 2 2dec2 x0 ~x þ ~xl þ Oððj~xj þ jljÞ3 Þ: ~x ! f~ ð~x; lÞ ¼ ~x þ l þ l þ 2 2 1e c 1  e2 c2 1  e2 c2

ð11Þ

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720 of~ ð0; 0Þ o~x

of~ ð0; 0Þ ol

o2 f~

705

2dec3 x

Since f~ ð0; 0Þ ¼ 0, ¼ 1, ¼ 1, o~x2 ð0; 0Þ ¼ 1e2 c20 6¼ 0, the fixed point ð~x; lÞ ¼ ð0; 0Þ is a fold bifurcation point for map (11). The number of fixed points is changed at ð~xðlÞ; ~y ðlÞÞ ¼ ð0; 0Þ as l ¼ 0, and if 1  c2 e2 < 0ð> 0Þ, then two new fixed points are created for l > 0 (l < 0) and disappear for l < 0(l > 0). From the above analysis, we have the proposition: Proposition 3. The map (2) undergoes a fold bifurcation at Z2 ðxðb0 Þ; yðb0 ÞÞ if the following conditions are satisfied sffiffiffiffiffiffiffiffiffiffiffi d 2ðc  1Þ 1 1 : c > 1; ec 6¼ 1; 2  dec þ 6¼ 0; and b0 ¼ ce 3 c Moreover, if 1  c2 e2 < 0 (resp. 1  c2 e2 > 0), then two fixed points bifurcate from Z2 for b < b0 (resp. b > b0 ), and coalesce as the fixed point Z2 at b ¼ b0 and disappear for b > b0 (resp. b < b0 ). We next consider the flip bifurcation of map (2). The characteristic equation associated with the linearization of the map (2) about the fixed point Zðx; yÞ is given by k2 þ pðbÞk þ qðbÞ ¼ 0;

ð12Þ

where pðbÞ ¼ 2 þ dec  dð1  x2 Þ=e; qðbÞ ¼ 1  dec þ dð1  x2 Þ=e  d2 cð1  x2 Þ þ d2 : Suppose b satisfies 1  pðb Þ þ qðb Þ ¼ 0:

ð13Þ

Correspondingly, the eigenvalues of the fixed point Zðxðb Þ; yðb ÞÞ are k1 ¼ 1 and k2 ¼ 1  pðb Þ. The condition jk2 j ¼ j1  pðb Þj 6¼ 1 leads to 2 þ dec  dð1  x2 ðb ÞÞ=e 6¼ 0; 2:

ð14Þ

~ ¼ b  b , we transform the fixed point Zðxðb Þ; yðb ÞÞ ¼ ðx0 ; y0 Þ of map (2) to Let u ¼ x  xðb Þ, v ¼ y  yðb Þ and b  ~ the origin, and take b as a new dependent variable, then map (2) becomes   1 0 0 1 u3 þ3u2 x0 þ3ux20 uþd u  v =e u 3 B C B ~ C C: ð15Þ @b A ! B ~ @ A b v  v þ deðu þ b~  cvÞ From (13) and (14), we have 2  dec 6¼ 0 and

ecðdec  2Þ  d 6¼ 0:

Provided 4  4dec þ d2 ðe2 c2  1Þ 6¼ 0; we construct an invertible matrix 0 1 1 1 d2 ðdec2Þ2 B 0 0 T ¼B d2 ðdec2Þ @ de dec2

eðdec2d2 Þ dð2decÞ

and use the translation 0 1 0 1 ~x u C ~ A ¼ T B @b l @ A; ~y v

eðdec2Þ d

ð16Þ 1 C C; A

706

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

the map (15) can be written 1 0 1 0 10 1 0 ~x ~x f ð~x; l; ~y Þ 1 1 0 A; @ l A ! @ 0 1 A@ l A þ @ 0 0 ~y ~y gð~x; l; ~y Þ 0 0 1  pðb Þ

ð17Þ

where f ð~x; l; ~y Þ ¼ 

ad ð3x0 ð~x þ l þ ~y Þ2 þ ð~x þ l þ ~y Þ3 Þ; 3e

gð~x; l; ~y Þ ¼ 

bd ð3x0 ð~x þ l þ ~y Þ2 þ ð~x þ l þ ~y Þ3 Þ; 3e



ðdec  2Þ2 ; 4  4dec þ d2 ðe2 c2  1Þ



d2 ðdec  2Þ ~ b: d  ðdec  2Þ2



d2 ; 4  4dec þ d2 ðe2 c2  1Þ

2

We again apply the center manifold theorem to determine the nature of the bifurcation of the fixed point ð~x; ~y Þ ¼ ð0; 0Þ at l ¼ 0. There exists a center manifold for (17) which can be represented as follows W c ð0Þ ¼ fð~x; l; ~y Þ 2 R3 j~y ¼ h ð~x; lÞ; h ð0; 0Þ ¼ 0; Dh ð0; 0Þ ¼ 0g: We assume h ð~x; lÞ ¼ b1~x2 þ b2~xl þ b3 l2 þ Oððj~xj þ jljÞ3 Þ:

ð18Þ

By the approximate computation for center manifold, we obtain     bdx0 2bdx0 1 bdx0 2 1 ; b2 ¼ þ 1 ; b3 ¼ þ þ1 : b1 ¼ pðb Þe pðb Þe pðb Þ pðb Þe p2 ðb Þ pðb Þ Thus, the map restricted to the center manifold is given by adx0 2 adx0 2 2adx0 ad ~x  ~xl  ð6b1 x0 þ 1Þ~x3 l  e e e 3e ad ad ad 2  ð2b1 x0 þ 2b2 x0 þ 1Þ~x l  ð2b2 x0 þ 2b3 x0 þ 1Þ~xl2  ð6b3 x0 þ 1Þl3 þ Oððj~xj þ jljÞ4 Þ: e e 3e

F : ~x ! ~x þ l 

In order for map (19) to undergo a flip bifurcation, we require that 0 ! !2 !1 2 3  oF o2 F o2 F  1 o F 1 o F A þ 6¼ 0 and a2 ¼ @ þ2 a1 ¼   2 2 3 ol o~x o~xol  2 o~x 3 o~x  ð0;0Þ It is obviously that a1 ¼

0  6adx e

ð19Þ

6¼ 0: ð0;0Þ

6¼ 0, and

2

a2 ¼

2ðdec  2Þ dðð3dc þ 5ecÞðdec  2Þ  3d2 ec þ dÞ : 3eð4  4dec þ d2 ðe2 c2  1ÞÞðecðdec  2Þ  dÞ

ð20Þ

If a21 ¼ ð3dc þ 5ecÞðdec  2Þ  3d2 ec þ d 6¼ 0;

ð21Þ

then a2 6¼ 0. Summarize the above analysis into the following proposition. Proposition 4. The map (2) undergoes a flip bifurcation at Zðxðb Þ; yðb ÞÞ if the conditions (13), (14), (16), and (21) are satisfied. Moreover, if a2 > 0 (resp. a2 < 0) in (20), then period-2 points that bifurcate from this fixed point are stable (resp. unstable). Finally, we give the condition of existence of Hopf bifurcation by using the Hopf bifurcation theorem in [8], here d is as bifurcation parameter.

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

707

The characteristic equation associated with the linearization of the map (2) at fixed point Zðxðd0 Þ; yðd0 ÞÞ ¼ ðx0 ; y0 Þ is given by k2 þ pðdÞk þ qðdÞ ¼ 0; where pðdÞ ¼ 2 þ dec  dð1  x20 Þ=e; qðdÞ ¼ 1  dec þ dð1  x20 Þ=e  d2 cð1  x20 Þ þ d2 : Correspondingly, the eigenvalues of the characteristic equation are pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pðdÞ  p2 ðdÞ  4qðdÞ ; k1;2 ¼ 2 where 2

2

"

p ðdÞ  4qðdÞ ¼ d

1  x20 þ ec e

2

# 4 :

The eigenvalues k1;2 are complex conjugate for p2 ðdÞ  4qðdÞ < 0, which leads to 1  eð2  ecÞ < x20 < 1 þ eð2 þ ecÞ:

ð22Þ

Let d0 ¼

ec  ð1  x20 Þ=e 1  cð1  x20 Þ

for 1  cð1  x20 Þ 6¼ 0;

ð23Þ

we have qðd0 Þ ¼ 1. We translate the fixed point Zðxðd0 Þ; yðd0 ÞÞ to the origin by the translation u ¼ x  x0 , v ¼ y  y0 . After the transformation, the map (2) becomes   !   u3 þ3u2 x0 þ3ux20 u  v =e uþd u 3 : ð24Þ G: ! v v þ deðu  cvÞ For d ¼ d0 , the eigenvalues of the matrix associated with the linearized map (2) at fixed point Oð0; 0Þ are complex conjugate with modulus 1 which are written as pðd0 Þ i pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  k; k ¼ 4qðd0 Þ  p2 ðd0 Þ; 2 2 where pðd0 Þ ¼ 2 þ d0 ec  d0 ð1  x20 Þ=e; qðd0 Þ ¼ 1  d0 ec þ d0 ð1  x20 Þ=e  d20 cð1  x20 Þ þ d20 : Under the conditions (22) and (23), there are jkj ¼ ðqðdÞÞ1=2

and



djkðdÞj 1 jd¼d0 ¼ d0 ð1  cð1  x20 ÞÞ 6¼ 0: dd 2

In addition, pðd0 Þ 6¼ 0; 1 leads to d0 ðec  ð1  x20 Þ=eÞ 6¼ 2; 3; n

then we have k ðd0 Þ 6¼ 1, n ¼ 1; 2; 3; 4. Let 0 ecð1x2 Þ=e 1 0 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 2 Þ=eÞ2 4ðecþð1x B C 0 T ¼@ A 2e pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 0 2 2 4ðecþð1x0 Þ=eÞ

ð25Þ

708

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

    x u , the map (24) becomes ¼T and use the translation y v pffiffiffiffiffiffiffiffiffiffiffiffiffi 1 0       i 4p2 ðd0 Þ pðd0 Þ x 2 2 A x þ f ðx; yÞ ; ffi ! @ pffiffiffiffiffiffiffiffiffiffiffiffi i 4p2 ðd0 Þ y y gðx; yÞ pðd0 Þ 2

ð26Þ

2

where f ðx; yÞ ¼ 0; 2

0

12

13 3

0

x20 Þ=e

x20 Þ=e

d6 B ec  ð1  1B ec  ð1  C C7 gðx; yÞ ¼ 4  @ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi x þ y A x0  @ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi x þ y A 5: e 3 2 2 2 2 4  ðec þ ð1  x0 Þ=eÞ 4  ðec þ ð1  x0 Þ=eÞ Notice that (26) is exactly in the form on the center manifold, in which the coefficient a [8] is given by " # 2 ð1  2kÞk 1 .11 .20  j.11 j2  j.02 j2 þ Re ðk.21 Þ; a ¼ Re 2 1k where 1 .20 ¼ ½fxx  fyy þ 2gxy þ iðgxx  gyy  2fxy Þ; 8 1 .11 ¼ ½ðfxx þ fyy Þ þ iðgxx þ gyy Þ; 4 1 .02 ¼ ½ðfxx  fyy  2gxy Þ þ iðgxx  gyy þ 2fxy Þ; 8 .21 ¼

1 ½ðfxxx þ fxyy þ gxxy þ gyyy Þ þ iðgxxx þ gxyy  fxxy  fyyy Þ: 16

Thus, an complex calculation gives "  d2 x20 ð1 þ s2 Þ 2dðec þ ð1  x20 Þ=eÞððpðd0 ÞÞ2 þ 1 þ 3pðd0 ÞÞ ððpðd0 ÞÞ2  3 þ pðd0 ÞÞð1  s2 Þ  a¼ 2 16e pðd0 Þ þ 2 #  dð1  d0 ecÞð1 þ s2 Þ jd¼d0 ;  3ð1 þ s2 Þ1  8e

ð27Þ

where dec  dð1  x20 Þ=e jd¼d0 : s ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 4qðd0 Þ  p2 ðd0 Þ From the above analysis, we have the Proposition 5. Proposition 5. The map (2) undergoes a Hopf bifurcation at fixed point Zðxðd0 Þ; yðd0 ÞÞ if the conditions (22), (23) and (25) hold and a 6¼ 0 in (27). Moreover, if a < 0 (resp. a > 0) and d ¼ djkðdÞj jd¼d0 > 0, then an attracting (resp. repelling) invariant dd closed curve bifurcates from the fixed point for d > d0 (resp. d < d0 ).

3. Existence of Marottos chaos In this section, we rigorously prove that the map (2) possesses chaotic behavior in the sense of Marotto’s definition [17].

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

709

We first present Marotto’s chaos definitions and theorem which are quoted from [17]. For any map F : Rn ! Rn , and any positive integer K, let F K represent the composition of F with itself K times. For a differentiable function F , let DF ðZÞ denote the Jacobian matrix of F evaluated at the point Z 2 Rn , and jDF ðZÞj its determinant. Let Br ðZÞ denote the closed ball in Rn of radius r centered at the point Z and B0r ðZÞ its interior. Also let kZk be the usual Euclidean norm of Z in Rn . Definition 1. Let F be differentiable in Br ðZ0 Þ. The point Z0 2 Rn is an expanding fixed point of F in Br ðZ0 Þ, if F ðZ0 Þ ¼ Z0 and all eigenvalues of DF ðZÞ exceed 1 in norm for all Z 2 Br ðZ0 Þ. Definition 2. Assume that Z0 is an expanding fixed point of F in Br ðZ0 Þ for some r > 0, then Z0 is said to be a snap-back  ¼ Z0 and jDF M ðZÞj  6¼ 0 for some positive integer M. repeller of F if there exists a point Z 2 Br ðZ0 Þ with Z 6¼ Z0 ; F M ðZÞ Theorem M. [17] If F possesses a snap-back repeller, then the map F is chaotic. That is , there exist (i) a positive integer N such that for each integer P P N , F has a point of period P ; (ii) a ‘‘scrambled set’’ of F , i.e, an uncountable set S containing no periodic points of F such that: (a) F ½S  S, (b) for every X ; Y 2 S with X 6¼ Y lim sup kF k ðX Þ  F k ðY Þk > 0

k!1

(c) for every X 2 S and any periodic point Y of F lim sup kF k ðX Þ  F k ðY Þk > 0

k!1

(iii) an uncountable subset S0 of S such that for every X ; Y 2 S0 : lim inf kF k ðX Þ  F k ðY Þk ¼ 0:

k!1

Now we theoretically give the condition of existence of chaotic phenomena for the map (2) in the sense of Marotto’s definition of chaos. Suppose Z0 ðx0 ; y0 Þ be the fixed point of map (2). We first give the conditions such that the Z0 is a snap-back repeller. The eigenvalues associated with the fixed point Z0 are given by pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pðx0 Þ  p2 ðx0 Þ  4qðx0 Þ ; k1;2 ¼ 2 where pðxÞ ¼ 2 þ dec  dð1  x2 Þ=e; qðxÞ ¼ 1  dce þ dð1  x2 Þ=e  d2 cð1  x2 Þ þ d2 : According to Definition 1, we begin to find a neighborhood Br ðZ0 Þ of Z0 in which the norms of conjugate complex eigenvalues exceed 1 for all Z 2 Br ðZ0 Þ. 2

Let s1 ðxÞ ¼ p2 ðxÞ  4qðxÞ ¼ de2 ½ðx2  1  e2 cÞ2  4e2 . It is easy to see that if 1  2e þ ce2 > 0; then the equation s1 ðxÞ ¼ 0 has four real roots pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi x1 ¼  1 þ 2e þ ce2 ; x2 ¼  1  2e þ ce2 ; pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi x3 ¼ 1  2e þ ce2 ; x4 ¼ 1 þ 2e þ ce2 : And s1 ðxÞ < 0 for x 2 A1 ¼ ðx1 ; x2 Þ [ ðx3 ; x4 Þ. Let s2 ðxÞ ¼ qðxÞ  1 ¼ de ½ðdec  1Þx2  dec þ 1  e2 c þ de.

ð28Þ

710

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

Under the conditions dec  1 6¼ 0

and 1 þ

eðec  dÞ > 0; dec  1

ð29Þ

the equation s2 ðxÞ ¼ 0 has two real roots sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi eðec  dÞ ; xþ ¼ 1 þ dec  1

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi eðec  dÞ x ¼  1 þ : dec  1

And s2 ðxÞ > 0 for all x 2 A2 ¼ ð1; x  [ ½xþ ; þ1Þ if dec  1 > 0, and x 2 A3 ¼ ½x ; xþ  if dec  1 < 0. Lemma 1 1e2 c 1 2 ; ec ; 1þec g < d < 1, and 1  2e þ ce2 > 0, then A1 \ A2 6¼ ;. (i) If c > 1, maxfðc1Þe 2

1e c 1 2 ; 0g < d < minfðc1Þe ; ec ; 1g, 1  2e þ ce2 > 0, or c < 1; 1  2e þ ce2 > 0; and0 < d < 1, then (ii) If c > 1, maxfec1 A1 \ A3 6¼ /. Moreover, if one of above conditions holds and the x-coordinate of fixed point Z0 ðx0 ; y0 Þ, x0 , satisfies x0 2 Ix0  A1 \ A2 ðA1 \ A3 Þ, then Z0 ðx0 ; y0 Þ is expanding fixed point of (2) in UZ0 ¼ fðx; yÞjx 2 Ix0 ; y 2 Rg.

Due to Definition 2 of snap-back repeller, we need to find one point Z 2 Br ðZ0 Þ such that Z 6¼ Z0 , F M ðZÞ ¼ Z0 , jDF M ðZÞj 6¼ 0 for some positive integer M. In fact, we have (

x þ dðx  x3 =3  yÞ=e ¼ x1 ; y þ deðx þ b  cyÞ ¼ y1

ð30Þ

x1 þ dðx1  x31 =3  y1 Þ=e ¼ x0 ; y1 þ deðx1 þ b  cy1 Þ ¼ y0 :

ð31Þ

and (

Now a F 2 map has been constructed to map the point Zðx; yÞ to the fixed point Z0 ðx0 ; y0 Þ after two iteration if there are solutions different from Z0 for Eqs. (30) and (31). By the calculation, the solutions different from Z0 for (31) are obtained as follows

(a)

(b)

Fig. 2. (a) Bifurcation diagram of fixed point of map (2) in (d–x) plane for b ¼ 0:6668, c ¼ 0:5, and e ¼ 0:01; (b) invariant cycle corresponding to d ¼ 0:052 in Fig. 2(a).

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

(a)

711

(c)

(b)

Fig. 3. Bifurcation diagrams of fixed points of map (2) for d ¼ 0:05, e ¼ 0:01 and different c: (a) c ¼ 0:5, (b) c ¼ 2 and (c) c ¼ 1:52.

beta 0.6710 0.6700 0.6690 0.6680 0.6670 0.6660 0.00

0.10

0.20

0.30

0.40

0.50 delta

Fig. 4. Hopf bifurcation locus for c ¼ 0:5, e ¼ 0:01.

Fig. 5. (a) Bifurcation diagram in (b–x) plane for c ¼ 2, d ¼ 0:05, e ¼ 0:01; (b) local amplification of (a) for showing the occurrences of period-6 and period-doubling reverse cascades into chaos, transient chaotic sets and interior crisis.

712

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

8 <

pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi de x0  3x20 þ12ð1þde þ1dec Þ

x1 ¼ 2 : y ¼ y  ðx1 x0 Þde 1 0 1dec

;

ð32Þ

de Þ. for x20 < 4ð1 þ de þ 1dec

Substituting (32) into (30) and solving x, we have !   e de 3b e de 3x0 3 x1 þ þ1þ xþ þ3  ¼ 0; x 3 d 1  dec c d ð1  decÞ2 cð1  decÞ2 where x1 ¼ x1 . We next need to find one real root x 2 Ix0 of (33). Let x ¼ s þ x0 , (33) becomes !   e de de e 3 2 2 þ 3ðx0  x1 Þ ¼ 0: s þ 3s x0 þ 3s x0   1   d 1  dec ð1  decÞ2 d

ð33Þ

ð34Þ

Fig. 6. Phase portraits of map (2) for c ¼ 2, d ¼ 0:05, e ¼ 0:01, and different b: (a) b ¼ 0:008, (b) b ¼ 0:01786, (c) b ¼ 0:018 and (d) b ¼ 0:025.

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

713

It is easy to see that if (34) has a real nonzero root s in a small neighborhood of zero, then x ¼ s þ x0 is just the 2 ð1decÞ2 required root of (33). Note that the constant item in (34), 3ðx0  x1 Þe d dð1decÞ 2 , may change from negative to positive 1 1 with respect to the parameter d. From d2  ð1  decÞ2 ¼ 0, we have d ¼ ceþ1 or d ¼ ce1 . Thus, there exists at least a d near d such that (33) has a real nonzero root x 2 Ix0 . de Lemma 2. If x20 < 4ð1 þ de þ 1dec Þ, then (33) has a real root x different from x0 , x 2 Ix0 for d near d .

From (30) and (31), we also obtain y ¼

y0 ð1  decÞ  ðx1  x0 Þde  deð1  decÞðx þ bÞ ð1  decÞ2 2

:

2

0Þ 0Þ Let UZ0 ¼ fðx; yÞj ðxx þ ðyy 6 1; jx  x0 j < rx0 ; x0  rx0 2 Ix0 ; ry0 ¼ jy   y0 j þ g; g is some positive constantg: r2 r2 x0

y0

Obviously, if the conditions in Lemmas 1 and 2 are satisfied and x0 2 Ix0 , then Z0 is a snap-back repeller in UZ0 . Thus, we have the theorem.

Fig. 7. (a) Bifurcation diagram in (d–x) plane for b ¼ 0:8, c ¼ 0:3, e ¼ 0:05; (b) maximum Lyapunov exponents corresponding to (a).

714

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

de Theorem A. If one of conditions in Lemma 1 holds, x0 2 Ix0 , and x20 < 4ð1 þ de þ 1dec Þ, then there exists at least a d near d such that Z0 ðx0 ; y0 Þ is a snap-back repeller of map (2), and hence map (2) is chaotic in the sense of Marotto’s definition.

We give specific values of the parameter for illustrating the condition in Lemmas 1 and 2 can be realized. Example. For c ¼ 0:8, b ¼ 0:45, e ¼ 0:05, d ¼ 0:96, there is the fixed point Z0 ðx0 ; y0 Þ ¼ ð0:9831063014; 0:6663828768Þ, and its eigenvalues are k ¼ 1:30242  i0:897464. From s1 ðxÞ ¼ 0 and s2 ðxÞ ¼ 0, we obtain x1 ¼ 1:04976, x2 ¼ 0:949737, x3 ¼ 0:949737, x4 ¼ 1:04976, and x ¼ 1:02364, thus A1 \ A3 6¼ ;. We take Br ðZ0 Þ ¼ fðx; yÞjðx  x0 Þ2 þ ðy  y0 Þ2 6 r2 ; r ¼ 0:02g. Let x1 ¼ x1 , (34) becomes f~ ðsÞ ¼ s3  3:306s  2:30005 ¼ 0. We can find one point Z ¼ ðx; yÞ where x ¼ 0:9832751197776846, y ¼ 0:6595259875424063 such that F 2 ðZÞ ¼ Z0 and DF 2 ðZÞ 6¼ 0. So, Z0 is a snap-back repeller.

Fig. 8. Phase portraits of map (2) for b ¼ 0:8, c ¼ 0:3, e ¼ 0:05, and different d: (a) d ¼ 0:051 and (b) d ¼ 0:0562.

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

715

4. Numerical simulations In this section, we present the results of the numerical simulations for the system (2). Let c ¼ 0:5, b ¼ 0:6668, e ¼ 0:01 and d is varying. Based on Proposition 5, after calculation for fixed point (x1 ; y 1 ) of map (2), we give ðx1 ; y 1 Þ ¼ ð1:00013; 0:666667Þ with D ¼ 5:0016 > 0, d ¼ 0:03166, a ¼ 3:76363, and d ¼ 0:01583. In fact, the fixed point (x1 ; y 1 ) is stable for d < 0:03166, and loses its stability at Hopf bifurcation value d ¼ 0:03166 which is labelled ‘‘1’’ in Fig. 2(a), as d increases, an invariant circle appears at d ¼ 0:052 which is shown in Fig. 2(b). Next b is taken as a bifurcation parameter. The bifurcation diagrams of fixed points of map (2) for d ¼ 0:05, e ¼ 0:01 and different c are given in Fig. 3(a)–(c) by using Auto97 [4]. In Fig. 3(a), we show that there are Hopf bifurcation (labelled ‘‘1’’) emerging from fixed point ðx1 ; y 1 Þ ¼ ð1:00023; 0:6667Þ at b1 ¼ 0:66689, and flip bifurcation (labelled ‘‘2’’) occurring at fixed point ðx2 ; y 2 Þ ¼ ð1:1833216; 0:631006Þ and b2 ¼ 0:867818, which check the correctness of Proposition 4. Fig. 3(b) shows that there are two fixed point curves for c ¼ 2, d ¼ 0:05, e ¼ 0:01, and there are three bifurcations: flip bifurcation (labelled ‘‘1’’) occurring at fixed point ðx1 ; y 1 Þ ¼ ð1:18332; 0:631006Þ as b1 ¼ 0:0786908, Hopf bifurcation (labelled ‘‘2’’) emerging from fixed point ðx2 ; y 2 Þ ¼ ð1:100015; 0:66667Þ as b2 ¼ 0:333183, and fold bifurcation (labelled ‘‘3’’) occurring at fixed point ðx3 ; y 3 Þ ¼ ð0:707107; 0:589256Þ and b3 ¼ 0:471405 which check the correctness of Propositions 3 and 4, respectively. When c ¼ 1:52, d ¼ 0:05 and e ¼ 0:01, the Hopf bifurcation and flip bifurcation separately occur on different fixed point curves which are shown in Fig. 3(c), a Hopf bifurcation (labelled

Fig. 9. (a) Bifurcation diagram pffiffi in (d–x) plane; (b) local amplification of (a) for d 2 ð0:7661; 0:7663Þ and (c) local amplification of (a) for d 2 ð0:796; 0:805Þ. Here b ¼ 33 þ 0:01, c ¼ 3, e ¼ 0:6.

716

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

‘‘1’’) occurs at ðx1 ; y 1 Þ ¼ ð1:00017; 0:66667Þ and b1 ¼ 0:013159 on the above one; and a flip bifurcation (labelled ‘‘3’’) is encountered at ðx3 ; y 3 Þ ¼ ð1:18332; 0:631Þ and b3 ¼ 0:224192 on the lower one. Fig. 4 is Hopf bifurcation locus in d–b space, we can observe that both d and b play important roles in Hopf bifurcation. Fig. 5(a) is the bifurcation diagram in (b–x) plane for e ¼ 0:01, c ¼ 2 and d ¼ 0:05. In fact, we can find that the fixed point ðx1 ; y 1 Þ ¼ ð1:18332; 0:631006Þ bifurcates period-two point at b1 ¼ 0:0786908, and inverse cascade of perioddoubling bifurcation occurs as b 2 ð0:031; 0:0786908Þ. When b 2 ð0:01; 0:031Þ, we show that period-window occurs in the chaotic region. The local amplification of Fig. 5(a) for b 2 ð0:0165; 0:0195Þ is plotted in Fig. 5(b) where we can see that there are period-6 orbits as b 2 ð0:0172; 0:01852Þ, interior crisis as b  0:0167, intermittency mechanic as b  0:01852, and transient chaos in period-window. The phase portraits at several critical values corresponding to Fig. 5 are shown in Fig. 6. The non-attracting chaotic sets at b ¼ 0:008, b ¼ 0:01786, and b ¼ 0:025 are shown as Fig. 6(a),(b) and (d), and corresponding Lyapunov exponents are f0:0027; 0:3137g, f0:0038; 0:2624g, and f0:0058; 0:2008g, respectively. The phase portraits of period-6 orbit at b ¼ 0:018 is plotted in Fig. 6(b).

Fig. 10. (a) pffiffi Phase portrait for period-26 (d ¼ 0:79728); (b) attracting invariant cycle (d ¼ 0:798) and (c) chaotic attractor (d ¼ 0:81). Here b ¼ 33 þ 0:01, c ¼ 3, e ¼ 0:6.

Fig. 11. Phase portraits of period-10 in two small periodic-windows at: (a) d ¼ 0:80159 and (b) d ¼ 0:80195. Here b ¼ e ¼ 0:6.

pffiffi 3 3

þ 0:01, c ¼ 3,

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

717

We next give the bifurcation diagram in (d–x) plane for c ¼ 0:3, e ¼ 0:05, b ¼ 0:8 in Fig. 7(a). The fixed point Z1 ðx1 ; y1 Þ ¼ ð1; 0:6667Þ of map (2) loses its stability at d ¼ 0:015 on account of the norm of its corresponding eigenvalues greater than 1, and there is an invariant circle when d > 0:015. The maximum Lyapunov exponents corresponding to bifurcation diagram Fig. 7(a) are shown in Fig. 7(b), where we can see that the maximum Lyapunov exponents are in the neighborhood of zero for d 2 ð0:012; 0:05Þ which is corresponding to quasiperiodic solutions or coexistence of chaos and quasiperiodic solutions. For d 2 ð0:05; 0:058Þ the maximum Lyapunov exponents are positive which corresponding to chaotic region. The phase portraits corresponding to d ¼ 0:051 and to d ¼ 0:0562 in the bifurcation diagram Fig. 7(a) are shown in Fig. 8(a) and (b), respectively. The Lyapunov exponents of the orbits in Fig. 8(a) and (b) are f0:0514; 0:660g and f0:552; 0:00g, respectively. Fig. 8(a) shows the coexistence of non-attracting chaotic set and periodic orbits, and Fig. 8(b) shows a non-attracting chaotic set. Finally we consider the cases in which there are three fixed points for pffiffi map (2) and d as the bifurcation parameter. Fig. 9(a) is the bifurcation diagram in (d–x) plane for c ¼ 3, b ¼ 33 þ 0:01 and e ¼ 0:6, and the initial values is ðx0 ; y0 Þ ¼ ð0; 0Þ. At this case, the map (2) has three fixed points: ðx1 ; y 1 Þ ¼ ð1:5429; 0:3185Þ; ðx2 ; y 2 Þ ¼ ð1:2346; 0:6073Þ and

ðx3 ; y 3 Þ ¼ ð0:3083; 0:2986Þ:

The bifurcation diagrams of local amplification of Fig. 9(a) are shown in Fig. 9(b) and (c) for 0:7661 6 d 6 0:7663 and 0:796 6 d 6 0:805, respectively. From Fig. 9(a), we observe that the orbit with initial value (0; 0) approaches to the stable fixed point (x1 ; y 1 ) for d < 0:743, and to the stable period-3 points for d 2 ð0:743; 0:755Þ, and there are two

Fig. 12. (a) Bifurcation diagram in (d–x) plane; (b) local amplification of (a) for d 2 ð0:8745; 0:878Þ; (c) local amplification of (a) for d 2 ð0:931; 0:935Þ and (d) local amplification of (a) for d 2 ð0:972; 0:974Þ. Here b ¼ 0:18, c ¼ 3, e ¼ 0:7.

718

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

different chaotic regions R1 ðd 2 ð0:764; 0:773ÞÞ and R2 ðd 2 ð0:79757; 0:823ÞÞ interspersed with several periodic windows. The chaotic region R1 with two coexisting chaotic attractors results from the period-doubling bifurcations in period-3, and disappears suddenly at d  0:773. The another chaotic attractor R2 with broaden size is from an invariant cycle bifurcation. In fact, at d ¼ 0:79757, a Hopf bifurcation occurs, and an attracting invariant cycle bifurcation from (x1 ; y 1 ) since a ¼ 0:40983 and d ¼ 2:05056 in Proposition 5. Fig. 10(b) presents the attracting cycle at d ¼ 0:798, the chaotic attractor at d ¼ 0:81 is shown in Fig. 10(c) and the associated Lyapunov exponents is computed as f0:336; 0:935g. The broadening chaotic attractor R2 disappears at d  0:823 when the new fixed point (x2 ; y 2 ) occurs (boundary crisis). The more complex bifurcation phenomena are exhibited in Fig. 9(b) and (c). From Fig. 9(b) it is clear that for d 2 ð0:76614; 0:76627Þ there are a period-15 orbits and some cascades of period-doubling bifurcation to chaos in the periodic window and the interior crisis at d  0:76627. From Fig. 9(c) we observe that there are the periodic-26 orbits at d ¼ 0:79728(the associated phase portrait is shown in Fig. 10(a)), period-10 orbits at d ¼ 0:80159 and at d ¼ 0:80195 in two small periodic windows, the associated phase portraits are given in Fig. 11(a) and (b), for d 2 ð0:8038; 0:80403Þ there are period-7 orbits and some cascades of doubling bifurcation to chaos in the periodic window, and the interior crisis occurs at d  0:80272 and d  0:80403.

Fig. 13. Phase portraits of map (2) for b ¼ 0:18, c ¼ 3, d ¼ 0:9322, e ¼ 0:7 and different initial values X0 : (a) X0 ¼ ð0; 0Þ and (b) X0 ¼ ð1; 0:6Þ.

Fig. 14. Phase portraits of map (2) for b ¼ 0:18, c ¼ 3, d ¼ 0:92 and e ¼ 0:7.

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

719

We also consider the case with three fixed points for another parameters c ¼ 3, b ¼ 0:18, and e ¼ 0:7. At this case, the map (2) has three fixed points ðx1 ; y 1 Þ ¼ ð1:4572; 0:4257Þ; ðx2 ; y 2 Þ ¼ ð1:3669; 0:5156Þ and

ðx3 ; y 3 Þ ¼ ð0:0904; 0:0901Þ:

The bifurcation diagram in (d–x) plane is shown in Fig. 12(a) for d 2 ð0:75; 1Þ, and the local amplificational bifurcation diagrams of Fig. 12(a) are shown in Fig. 12(b), (c) and (d) for d 2 ð0:8745; 0:878Þ, d 2 ð0:931; 0:935Þ and d 2 ð0:972; 0:974Þ, respectively. From Fig. 12(a) we observe that there are two different chaotic region R1 ðd 2 ð0:87; 0:9334ÞÞ and R2 ðd 2 ð0:972; 1ÞÞ, and several periodic windows in the chaotic regions. The two chaotic regions with interior crisis result from the period-doubling bifurcations in period-3, and there is a region for the coexistence of stable period-3 orbits with chaotic attractor (see Fig. 12(c) for d 2 ð0:9317; 0:93345Þ),the phase portraits with chaotic attractor and period-3 orbits are separately given in Fig. 13(a) and (b). Moreover, we find the Marotto’s chaotic attractor rests in the chaotic region R1 . In fact, when d ¼ 0:92, fixed point  1 ; y 1 Þ of map (2) have complex conjugate eigenvalues with norm greater than 1 that are kZ ¼ 0:7043  i0:8913. Zðx 2 2 From the analysis in Section 3, the elliptic region of Z is U  ¼ fðx; yÞj ðxx21 Þ þ ðyy21 Þ 6 1; rx ¼ 0:02; ry ¼ 0:7g in which Z Z

rx

ry

is expanding, and there exists a point Z  ¼ ðx ; y  Þ ¼ ð1:47374; 1:11497Þ such that F 2 ðZ  Þ ¼ Z and j DF 2 ðZ  Þ j6¼ 0. Thus Z is a snap-back repeller. The Marrotto’s chaotic attractor is shown in Fig. 14 for d ¼ 0:92, while the associated Lyapunov exponents are f0:34615; 0:22876g. On the other hand, the coexisting two chaotic attractors and several period-windows are shown in Fig. 12 (b) and (d), which are period-15 and period-9 windows, and associated perioddoubling cascades for 0:875 < d < 0:8753 and 0:8755 < d < 0:8775, respectively (see Fig. 12(b)), period-15 and 21 windows and associated period-doubling cascades for 0:9726 < d < 0:9732 and 0:9736 < d < 0:9738, respectively (see Fig. 12(d)).

References [1] Beltrami E. Mathematics for dynamic modeling. Springer-Verleg; 1997. [2] Braaksma B, Grasman J. Critical dynamics of the Bonhoeffer–Van der Pol equation and its chaotic response to periodic stimulation. Physica D 1993;68:265–80. [3] Brown D, Foweraker JPA, Marrs RW. Dynamic equilibria and oscillations of a periodically stimulated excitable system. Chaos, Solitions & Fractals 1995;5(3/4):359–69. [4] Doedel EJ, Fairgrient FT, Wang X. AUTO 97 Continuation and bifurcation software for ordinary differential equation, 1997. [5] FitzHugh R. Impulses and physiological states in theoretical models of nerve membrane. Biophys J 1961;1:445–66. [6] Flores G. Stability analysis for the slow traveling pulse of the FitzHugh–Nagumo system. SIAM J Math Anal 1991;22(2):392–9. [7] Freitas P, Rocha C. Lyapunov functional and stability for FitzHugh–Nagumo systems. J Different Equat 2001;169:208–27. [8] Guckenheimer J, Holmes. Nonlinaer oscillations, dynamical systems and bifurcations of vector fields. New York: Springer-Verlag; 1983. [9] Hale JK, Kocak H. Dynamics and bifurcations. Springer-Verlag; 1991. [10] Hassard BD. Bifurcation of periodic solutions of the Hodgkin–Huxley model for the squid giant axon. J Theoret Biol 1978;71:401–20. [11] Hayashi M. A note on the uniqueness of the closed orbit of the FithHugh–Nagumo system. Quart Appl Math 2000;LViii(1):171–6. [12] Hayashi M. Non-existence of homoclinic orbits and global asymptotic stability of FitzHugh–Nagumo system. Vietnam J Math 1999;27(4):335–43. [13] Hayashi M. On the uniqueness of the closed orbit of FitzHugh–Nagumo system. Math Japonica 1997;46(2):331–6. [14] Hodgkin AL, Huxley AF. A quantitative description of membrane current and its application to conduction and excitation in nerve. J Physiol 1952;117:500–44. [15] Jing ZJ, Jia ZY, Wang RQ. Chaos behavior in the discrete BVP oscillator. Int J Bifurcat Chaos 2002;12(3):619–27. [16] Jones CKRT. Stability of travelling wave solution of the FitzHugh–Nagumo system. Trans Am Math Soc 1984;286(2):431–69. [17] Marotto FR. Snap-back repellers imply chaos in Rn . J Math Anal Appl 1978;63:199–223. [18] Mattleig RMM, Molenaar. Ordinary differential equations in theory and practice. John Wiley & Sons Ltd; 1996. [19] Mira C, Gardini L, Barugola A, Cathala J. Chatic dynamics in two-dimensional noninvertible maps. Singapore: World Scientific; 1996. [20] Murray JD. Mathematical biology. Springer-Verlag; 1993. [21] Nagumo J, Arimoto S, Yoshizawa S. An active pulse transmission line simulating nerve axon. Proc IRE 1962;50:2061–70. [22] Ranch J, Smoller J. Qualitative theory of the FitzHugh–Nagumo equations. Adv Math 1978;27:12–44. [23] Rinzel J, Miller RN. Numerical calculation of stable and unstable periodic solutions to the Hodgkin–Huxley equations. Math Biosci 1980;49:27–59. [24] Rocsoreanu C, Giurgiteanu N, Georgescu A. Connections between saddles for the FitzHugh–Nagumo system. Int J Bifurcat Chaos 2001;11(2):533–40.

720

Z. Jing et al. / Chaos, Solitons and Fractals 21 (2004) 701–720

[25] Seydel R. Practical bifurcation and stability analysis-from equilibrium to chaos. Springer-Verlag; 1994. [26] Smoller J. Shock waves and reaction diffusion equations. Springer-Verlag; 1994. [27] Van der Pol B, Van der Mark J. The heart beat considered as a relaxation oscillation and an electrical model of the heart. Philos 1928;6(May):763–75. [28] Wiggins S. Introduction to applied nonlinear dynamical systems and chaos. Springer-Verlag; 1990. [29] Winfree AT. Varieties of spiral wave behavior: An experimentalist’s approach to the theory of excitable media. Chaos 1991;1:303– 34. [30] Zhang H, Holdem AV. Chaotic meander of spiral waves in the FitzHugh–Nagumo system. Chaos, Solitons & Fractals 1995;5:661–70.