3 June 2002
Physics Letters A 298 (2002) 109–116 www.elsevier.com/locate/pla
Chaotic attractors in delayed neural networks Hongtao Lu Department of Computer Science and Engineering, Shanghai Jiaotong University, Shanghai, 200030, PR China Received 4 July 2001; received in revised form 17 April 2002; accepted 22 April 2002 Communicated by A.R. Bishop
Abstract This Letter investigates the complex dynamical behavior of delayed neural networks with two neurons with the help of computer simulations. It has been shown that such networks may exhibit chaotic dynamics undergoing a period-doubling bifurcation process. In some parameter domains, interesting phenomena of coexistence of periodic orbits and chaotic attractors have been observed. 2002 Elsevier Science B.V. All rights reserved. Keywords: Chaos; Neural networks; Delay; Bifurcation
1. Introduction Hopfield realized that in hardware implementation of neural network models, time delays occur due to finite switching speed of the amplifiers [1]. Marcus and Westervelt [2] first proposed a neural network model with delay. Another type of delayed neural networks is delayed cellular neural networks (DCNN’s) introduced by Roska and Chua [3]. There has recently been increasing interest in the dynamics of neural networks with delays. Most of existing works focused on the stability analysis and periodic oscillation of this kind of neural networks [2,4–8]. Stability and complex dynamical behaviour of DCNN’s have been investigated [9–15]. Gilli has reported some strange attractors in DCNN with two neurons [10]. In this Letter we will investigate complex dynamics of delayed neural networks with the help of numerical simulations. We involve only networks with two neurons and single deE-mail address:
[email protected] (H. Lu).
lay, variety of dynamics have been found, such as periodic oscillations with different periods, quasi-periodic solutions, and chaotic attractors. Some interesting phenomena of coexistence of periodic orbits and chaotic attractors have been observed.
2. The network description The neural networks with single delay considered in this Letter are described by the following differential equations with delay: dxi (t) = −ai xi (t) + bij f xj (t) dt n
j =1
+
n
cij f xj (t − τ )
j =1
+ Ii ,
i = 1, 2, . . . , n,
(1)
where all ai > 0, bij and cij ’s are real numbers, τ 0 represents the time delay, Ii external inputs.
0375-9601/02/$ – see front matter 2002 Elsevier Science B.V. All rights reserved. PII: S 0 3 7 5 - 9 6 0 1 ( 0 2 ) 0 0 5 3 8 - 8
110
H. Lu / Physics Letters A 298 (2002) 109–116
The input–output transfer function f (x) is chosen as f (x) = tanh(x). We consider a two-neuron network, in this case the model can be rewritten as dx(t) = −Ax(t) + Bf x(t) + Cf x(t − τ ) + I, dt (2) where
T x(t) = x1 (t), x2 (t) , T f x(t) = f x1 (t) , f x2 (t) , T f x(t − τ ) = f x1 (t − τ ) , f x2 (t − τ ) , I = (I1 , I2 )T , and
A= C=
a1 a2 c11 c21
,
B=
b11 b21
b12 , b22
c12 . c22
The difference differential equation (2) can be categorized as a kind of functional differential equations [16], and rewritten as dx(t) = −Axt (0) + Bf xt (0) + Cf xt (−τ ) + I, dt (3) where xt is a continuous mapping defined on [−τ, 0] as xt (θ ) = x(t + θ ), the right-hand side of Eq. (3) defines a functional mapping C([−τ, 0], R 2 ) to R 2 , where C([−τ, 0], R 2 ) denotes the set of all continuous mappings from [−τ, 0] to R 2 . The solution space of Eq. (3) is infinite-dimensional, with initial conditions as any continuous functions defined on the closed interval [−τ, 0].
3. Equilibrium points and eigenvalues An equilibrium point of Eq. (3) is a vector xe = (xe1 , xe2 )T that satisfies condition −Axe + (B + C)f (xe ) + I = 0. In the remainder of this Letter, we always assume I = 0. Since the function f (x) = tanh(x) is odd symmetry, it is easy to see that the origin (0, 0) is always an equilibrium point, and if xe is an equilibrium point of (3), so is −xe . This means equilibrium points of (3)
are distributed in the real plane symmetrically with respect to the origin (0, 0). Qualitative behavior of Eq. (3) depends on the behavior near its equilibrium points, and solutions of a functional differential equation have the same qualitative behavior near an equilibrium point as the solutions of its linearized equation at this equilibrium, if no eigenvalues of the linearized equation lie on the imaginary axis (see [16, Chapter 10]). The behavior of a linear equation is determined essentially by the eigenvalues and corresponding eigenvectors of its characteristic equation. Thus, we first consider the linearized equation of (3) at xe , which is given as dx(t) = −Ax(t) + BDx(t) + CDx(t − τ ), dt where 1 − (tanh(xe1 ))2 D=
(4)
1 − (tanh(xe2 ))2
.
Eigenvalues are the roots of the transcendental equation: λE + A − BD − CDe−λτ = 0, where 1 E=
(5)
1
.
An eigenvector v = (v1 , v2 )T corresponding to a given eigenvalue λ0 is a solution of the linear equation: λ0 E + A − BD − CDe−λ0 τ v = 0.
(6)
Usually, Eq. (5) has infinite number of roots, there exists no analytic method to solve this type of equations, in the following we will find its roots numerically. If real parts of all roots are negative, then the equilibrium point xe is local stable [16], i.e., any solution started nearby will approach to xe . If there exist roots with positive real parts, the equilibrium xe is local unstable, trajectories will converge to xe along the direction of the eigenvectors corresponding to eigenvalues with negative real parts, and diverge from xe along the direction of the eigenvectors belong to eigenvalues with positive real parts.
H. Lu / Physics Letters A 298 (2002) 109–116
4. Bifurcation and chaos In what follows we fixed some parameters of Eq. (2) as: 1 A= , 1 the delay τ = 1 and I = (0, 0)T , and use the fourth order Runge–Kutta method to solve Eq. (2), the terms with delay are treated as input to the system. The stepsize is chosen as 0.005, and without loss of generality, initial conditions are always chosen as constant functions on [−τ, 0]. We first choose parameters as −1.5 −0.1 2.0 −0.1 , C= , B= −0.2 c22 −5.0 b22 and only b22 and c22 are allowed to change. Some dynamical behavior are shown in Figs. 1, 2 and 3(a), where parameters b22 and c22 are varied, while keeping a relation b22 + c22 = 0.5, this makes the system have identical equilibrium points in all figures. There are total three equilibrium points calculated numerically as x−1 = (−0.37614, 2.35978), x0 = (0, 0), x1 = (0.37614, −2.35978), and marked by bullets in these figures. Fig. 1(a) shows period 1 solution when b22 = 0.3, c22 = 0.2, in this case there exist two separate periodic orbits when started from different initial conditions. One periodic orbit oscillates around x1 , and another around x−1 . As increasing parameter b22 and decreasing c22 , solutions with period 2 and 4 appear, Fig. 1(b) and (c) show a pair of period two and four orbits with parameters as b22 = 0.5, c22 = 0 and b22 = 0.8, c22 = −0.3, respectively. As further change of parameters, the system enters chaotic domain by a route of period-doubling bifurcation process. Fig. 1(d) shows two coexisting separate chaotic attractors, oscillating around x−1 and x1 , with different initial conditions and parameters as b22 = 1.5, c22 = −1.0. Note that there are two separate single-scroll-like [17] attractors in the figures, although they look like a whole one. As the parameters change slightly, two separate single-scroll-like attractors merge to a double-scrolllike [17] one oscillating around x1 and x−1 alternatively, such an attractor is displayed in Fig. 2(a) with b22 = 2.0, c22 = −1.5. Fig. 2(b) depicts a fully developed double-scroll-like chaotic attractor with b22 =
111
3.0, c22 = −2.5. With further increasing of b22 and decreasing of c22 , chaotic and periodic solutions appear alternatively, Figs. 2(c), (d) and 3(a) illustrate a period orbit when b22 = 4.0, c22 = −3.5, a chaotic attractor with b22 = 4.5, c22 = −3.5, and again a periodic orbit with b22 = 5.0, c22 = −4.5. From these figures, we can see that as parameters vary, the system enters chaos undergoing a period-doubling process, and dynamics of the system is symmetry with respect to the origin due to the odd symmetry property of the vector field. One should also note that in all 9 figures, the equilibrium points of the network are fixed due to the special choice of the parameters, only eigenvalues and eigenvectors are changed, which are calculated numerically according to Eqs. (5) and (6), and listed in Table 1. Since the odd symmetry of the vector field, it is easy to prove that x1 and x−1 have identical eigenvalues and eigenvectors. From Table 1 we can see that for all cases all three equilibrium points are saddle with eigenvalues with both positive and negative real parts. Nonzero imaginary parts of eigenvalues at x1 and x−1 make it possible to exhibit complex behavior. In the following, we investigate the system behavior, fixing parameters as 2.0 −0.1 −1.5 −0.1 B= , , C= −5.0 4.5 −0.2 c22 and change only c22 . When taking c22 = −2.0, three coexisting period 1 orbits oscillating around three equilibrium points x−1 = (−0.3821, 3.8939), x0 = (0, 0), x1 = −x−1 are found. The periodic orbits have their own basins of attraction which do not overlap with each other, Fig. 3(b) shows such a case. When c22 = −2.5, there is only one global periodic orbit around x0 left, two other ones disappear, any initial condition will converge to this periodic orbit, no other attracting sets exist, Fig. 3(c) shows the case. Figs. 3(d) and 4(a) depict an interesting phenomenon of coexistence of two types of chaotic attractors at a same set of parameters when c22 = −3.0. Fig. 3(d) is one type of chaotic attractor only oscillating around x0 in a small area with simpler structure, while Fig. 4(a) another type oscillating all three equilibrium points with more complicated structure. There exists another type of coexistence of attractors: two small chaotic attractors and a large periodic orbit shown in Fig. 4(b) when c22 = −4.45, the chaotic attractors oscillate around x−1 = (−0.3659, 1.8699) and x1 = −x−1 in
112
H. Lu / Physics Letters A 298 (2002) 109–116
Fig. 1. (a) Period one orbits with parameters b22 = 0.2, c22 = 0.3 and initial conditions x1 (t) ≡ 0.4, x2 (t) ≡ 0.6 and x1 (t) ≡ −0.4, x2 (t) ≡ −0.6 for t ∈ [−1, 0], respectively. (b) Period two orbits with parameters b22 = 0.5, c22 = 0 and initial conditions x(t) ≡ 0.4, y(t) ≡ 0.6 and x(t) ≡ −0.4, y(t) ≡ −0.6 for t ∈ [−1, 0], respectively. (c) Period four orbits with parameters b22 = 0.8, c22 = −0.3 and initial conditions x(t) ≡ 0.4, y(t) ≡ 0.6 and x(t) ≡ −0.4, y(t) ≡ −0.6 for t ∈ [−1, 0], respectively. (d) Two separate single-scroll-like chaotic attractors with parameters b22 = 1.5, c22 = −1.0 and initial conditions x(t) ≡ 0.4, y(t) ≡ 0.6 and x(t) ≡ −0.4, y(t) ≡ −0.6 for t ∈ [−1, 0], respectively.
Fig. 2. (a) Double-scroll-like chaotic attractor with parameters b22 = 2.0, c22 = −1.5 and initial conditions x(t) ≡ 0.4, y(t) ≡ 0.6. (b) Fully developed double-scroll-like chaotic attractors with parameters b22 = 3.0, c22 = −2.5 and initial conditions x(t) ≡ 0.4, y(t) ≡ 0.6. (c) Periodic orbit with parameters b22 = 4.0, c22 = −3.5. (d) Chaotic attractors with parameters b22 = 4.5, c22 = −4.0.
H. Lu / Physics Letters A 298 (2002) 109–116
113
Fig. 3. (a) Periodic orbit with parameters b22 = 5.0, c22 = −4.5. (b) Three coexisting periodic orbits with parameters b22 = 4.5, c22 = −2.0 and started from initial conditions x(t) ≡ 0.4, y(t) ≡ 0.6, x(t) ≡ −0.4, y(t) ≡ −0.6 and x(t) ≡ 2.0, y(t) ≡ −3.0, respectively. (c) Periodic one orbit oscillating around (0,0) with parameters b22 = 4.5, c22 = −2.5. (d) One of two different coexisting chaotic attractors with parameters b22 = 4.5, c22 = −3.0 and initial condition x(t) ≡ 0.4, y(t) ≡ 0.4.
Fig. 4. (a) One of two different coexisting chaotic attractors with parameters b22 = 4.5, c22 = −3.0 and initial condition x(t) ≡ −0.5, y(t) ≡ 3.0. (b) Three types of coexisting attractors at parameters b22 = 4.5, c22 = −4.45, started from initial conditions x(t) ≡ 0.4, y(t) ≡ 0.6 (the periodic orbit) and x(t) ≡ −0.4, y(t) ≡ 3.0 or x(t) ≡ 0.4, y(t) ≡ −3.0 (the two chaotic attractors). (c) Enlargement of the chaotic attractor based at x(t) ≡ −0.4, y(t) ≡ 3.0. (d) Basins of attraction of three types of attractors in (b) and (c).
114
H. Lu / Physics Letters A 298 (2002) 109–116
Table 1 Dominant eigenvalues and corresponding eigenvectors of the linearized system (4) at x1 , x0 and x−1 x−1 , x1 Eigenvalues Fig. 1(a) Fig. 1(b) Fig. 1(c) Fig. 1(d) Fig. 2(a) Fig. 2(b)
Fig. 2(c) Fig. 2(d)
Fig. 3(a)
−0.9301 0.08264 ± 0.98511i −0.9466 0.08257 ± 0.98502i −0.96401 0.08249 ± 0.98492i −1.00896 0.08233 ± 0.9847i −1.0456 0.0822 ± 0.9845i −3.341 −1.1353 0.08195 ± 0.9842i −1.26363 0.08168 ± 0.9838i −2.4605 −1.3591 0.08155 ± 0.9837i −2.1351 −1.5155 0.08142 ± 0.9835i
x0
Eigenvectors
Eigenvalues
Eigenvectors
(1, −132.32) (1, −2.16827 ± 2.13692i) (1, −133.84) (1, −2.18589 ± 2.14711i) (1, −135.45) (1, −2.20377 ± 2.15737i) (1, −139.67) (1, −2.24657 ± 2.1816i) (1, −143.15) (1, −2.2781 ± 2.1991i) (1, −320.14) (1, −151.83) (1, −2.34366 ± 2.2347i) (1, −164.51) (1, −2.4128 ± 2.271i) (1, −271.51) (1, −174.02) (1, −2.44875 ± 2.2895i) (1, −246.52) (1, −189.57) (1, −2.48572 ± 2.308i)
0.8337 −0.2737 ± 1.0492i 0.901
(1, −3.3836) (1, −0.6726 ± 3.5289i) (1, −3.6276)
0.9769 −2.1721 1.205 −0.901 1.4327 −0.6704 2.1325 −0.4599
(1, −3.9348) (1, −10.2212) (1, −5.036) (1, −5.1764) (1, −6.3837) (1, −4.2712) (1, −11.714) (1, −3.5449)
3.0831 −0.3547 3.5849 −0.3188
(1, −20.575) (1, −3.2317) (1, −25.556) (1, −3.1338)
4.0879 −0.2897
(1, −30.617) (1, −3.0579)
Fig. 5. (a) Three periodic orbits coexist at parameters b22 = 4.5, c22 = −4.6. (b) Periodic one orbit at parameters b22 = 4.5, c22 = −6.0. −1.5 0.1 (c) Attractor at parameters b22 = 4.5, c22 = −8.5. (d) Similar chaotic attractor to Fig. 2(b) at parameters B = 2.0 0.1 , C = . 5.0 3.0
0.2 −2.5
H. Lu / Physics Letters A 298 (2002) 109–116
115
Table 2 Eigenvalues and eigenvectors of the linearized system (4) at x1 , x0 and x−1 for Figs. 5(d), 6(a),(b) x−1 , x1 Eigenvalues Fig. 5(d)
Fig. 6(a)
Fig. 6(b)
x0 Eigenvectors
Eigenvalues
x1 = (0.37614, 2.35978), x−1 = −x1 −3.341 (1, 320.14) −1.1353 (1, 151.83) 0.08195 ± 0.9842i (1, 2.34366 ∓ 2.2347i) x1 = (2.35978, 0.37614), x−1 = −x1 −3.341 (320.14, 1) −1.1353 (151.83, 1) 0.08195 ± 0.9842i (2.34366 ∓ 2.2347i, 1) x1 = (2.35978, −0.37614), x−1 = −x1 −3.341 (−320.14, 1) −1.1353 (−151.83, 1) 0.08195 ± 0.9842i (−2.34366 ± 2.2347i, 1)
Fig. 6. (a) Chaotic attractor at parameters B = −2.5 −0.2 C= .
3.0 5.0 , C = 0.1 2.0
Eigenvectors
2.1325 −0.4599
(1, 11.714) (1, 3.5449)
2.1325 −0.4599
(11.714, 1) (3.5449, 1)
2.1325 −0.4599
(−11.714, 1) (−3.5449, 1)
−2.5 0.2 . (b) Chaotic attractor at parameters B = 0.1 −1.5
3.0 −5.0 , −0.1 2.0
−0.1 −1.5
a small area, respectively, while the periodic orbit around x0 in a large area. Fig. 4(c) is the magnification of one of the chaotic attractors. Fig. 4(d) displays the basins of attraction of different attractors when initial conditions confined to constant functions on [−τ, 0], taking values in −1.5 x1 1.5 and −4.0 x2 4.0. Solutions started from the area marked by stars will converge to the small chaotic attractor around x−1 , those based at points in the area marked by circles will be attracted to the neighborhood of x1 , the remained area is the basin of attraction of the periodic orbit in Fig. 4(b) and (c). When c becomes −4.6, two chaotic attractors turn to period 1 orbits, and as c decreases further to −6.0, the two small periodic orbits disappear and only a large periodic orbit left as the global attracting limit set of the system, these are displayed in Fig. 5(a) and (b). When c reduces to
−8.5, the network has a quasi-periodic attractor shown in Fig. 5(c). The dynamics of the network have some symmetrical properties, in order to illustrate these properties, Figs. 5(d), 6(a) and (b) give three chaotic attractors with three sets of parameters as 2.0 0.1 −1.5 0.1 B= , C= , 5.0 3.0 0.2 −2.5 3.0 5.0 −2.5 0.2 B= , C= , 0.1 2.0 0.1 −1.5 and B=
3.0 −5.0 , −0.1 2.0
C=
−2.5 −0.2 , −0.1 −1.5
respectively. They can be compared with Fig. 2(b). Their equilibrium points, eigenvalues and eigenvectors
116
H. Lu / Physics Letters A 298 (2002) 109–116
are listed in Table 2. It can be seen that equilibrium points of Figs. 2(b), 5(d) and 6(a), (b) are located in the four quadrants of the x1 –x2 plane symmetrically, and all with identical eigenvalues, while eigen-directions (eigenvectors) are also symmetric.
tion Commission. Basic is affiliated to the Department of Computer Science of Shanghai Jiaotong University. We wish to thank the referees for their valuable suggestions and comments.
References 5. Conclusions Complex dynamical behavior of delayed neural network with two neurons have been investigated with help of numerical simulation. A period-doubling route to chaos has been observed and some interesting phenomena of coexistence of periodic orbits and chaotic attractors have been found. It is suggested that introduction of delays into neural network models makes their behavior more complicated. Here we fix all delays to 1, and just change some parameters, complex dynamics under changes of delays and rigorous treatment of bifurcation and chaos deserve more further investigation.
Acknowledgements This work is supported by the Natural Science Foundation of China under grant No. 69901002 and supported in part by BASICS, Center of Basic Studies in Computing Science, sponsored by Shanghai Educa-
[1] J.J. Hopfield, Proc. Nat. Acad. Sci. U.S.A. 81 (1984) 3088. [2] C.M. Marcus, R.M. Westervelt, Phys. Rev. A 39 (1989) 229. [3] T. Roska, L.O. Chua, Int. J. Circuit Theory Appl. 20 (1992) 469. [4] K. Gopalsamy, X.Z. He, Physica D 76 (1994) 344. [5] J. Belair, S.A. Campbell, P. van den Driessche, SIAM J. Appl. Math. 56 (1996) 245. [6] P. Baldi, A.F. Atiya, IEEE Trans. Neural Networks 5 (1994) 612. [7] Y.J. Cao, Q.H. Wu, IEEE Trans. Neural Networks 7 (1996) 1533. [8] Y. Zhang, Int. J. Systems Sci. 27 (2) (1996) 227. [9] P.P. Civalleri, M. Gilli, L. Pandolfi, IEEE Trans. Circuit Systems I 40 (1993) 157. [10] M. Gilli, IEEE Trans. Circuit Systems I 40 (1993) 849. [11] M. Gilli, IEEE Trans. Circuit Systems I 41 (1994) 518. [12] J. Cao, D. Zhou, Neural Networks 11 (1998) 1601. [13] J. Cao, Phys. Lett. A 261 (1999) 303. [14] J. Cao, Phys. Lett. A 267 (2000) 312. [15] J. Cao, Phys. Lett. A 270 (1999) 157. [16] J. Hale, Theory of Functional Differential Equations, Springer, New York, 1977. [17] L.O. Chua, M. Komuro, T. Matsumoto, IEEE Trans. Circuit Systems 33 (1986) 1073.