- Email: [email protected]
Effect of random noise on chaotic motion of a particle in a ϕ6 potential
Chaos, Solitons and Fractals 27 (2006) 127–138 www.elsevier.com/locate/chaos
Effect of random noise on chaotic motion of a particle in a /6 potential Zhongkui Sun a, Wei Xu a
a,*
, Xiaoli Yang
a,b
Department of Applied Mathematics, Northwestern Polytechnic University, Xi’an 710072, PR China b Department of Mathematics, Shaan’xi Normal University, Xi’an 710062, PR China Accepted 1 February 2005
Abstract The chaotic behaviors of a particle in a triple well /6 potential possessing both homoclinic and heteroclinic orbits under harmonic and Gaussian white noise excitations are discussed in detail. Following Melnikov theory, conditions for the existence of transverse intersection on the surface of homoclinic or heteroclinic orbits for triple potential well case are derived, which are complemented by the numerical simulations from which we show the bifurcation surfaces and the fractality of the basins of attraction. The results reveal that the threshold amplitude of harmonic excitation for onset of chaos will move downwards as the noise intensity increases, which is further verified by the top Lyapunov exponents of the original system. Thus the larger the noise intensity results in the more possible chaotic domain in parameter space. The effect of noise on Poincare maps is also investigated. Ó 2005 Published by Elsevier Ltd.
1. Introduction Collective investigations have been dedicated to the dynamics of a harmonically excited particle in a /4 potential [1–4] 1 1 V ðxÞ ¼ bx2 þ cx4 ; 2 4 where b and c are constants. Nevertheless, with the advent of the study of chaotic motion by means of strange attractors, top Lyapunov exponents, Poincare maps and fractal basin boundaries, it has become necessary to seek for better understanding of nonlinear system with higher order nonlinear terms. In this paper, we consider the extended Duffing oscillator by introducing a quintic term described by the following equation: x_ ¼ y; ð1Þ y_ ¼ bx cx3 dx5 ay þ c cos Xt þ rnðtÞ;
*
Corresponding author. Tel.: +86 29 88495453; fax: +86 29 88494174. E-mail addresses: [email protected] (Z. Sun), [email protected] (W. Xu).
0960-0779/$ - see front matter Ó 2005 Published by Elsevier Ltd. doi:10.1016/j.chaos.2005.02.033
128
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
Fig. 1. Phase portrait and potential of system (2) with b = 1.0, c = 0.5, d = 0.038. (a) Phase portrait and (b) potential.
where b, c and d are constants parameters; a is the damping ratio; n(t) represents a standard Gaussian white noise with Dirac-Delta function as its correction, namely, En(t) = 0, En(t)n(t + s) = d(s) and c and r denote the amplitude or level of the harmonic or noise excitation, respectively. The interest in such a mode is due to the fact that it is a universal nonlinear differential equation, and many nonlinear oscillators in physical, engineering and biological problems can be really described by the model or analogous ones. The potential from which is a /6 potential given by 1 1 1 V ðxÞ ¼ bx2 þ cx4 þ dx6 . 2 4 6 Depending on the set of the parameters, it can be considered at least three physically interesting situations where the potential is (i) single well, (ii) double well, and (iii) triple well. Throughout this paper, our analysis is carried out on the case of triple well with a double hump which does not lead to unbounded motion (see Fig. 1(b)). The effects of excitation, consisting of harmonic and random noise excitations, on nonlinear dynamical systems exhibiting chaotic behavior have been examined in recent years [5–12]. Especially the chaotic behaviors of this type of nonlinear system under periodic excitation have been extensively investigated (for instance, see Refs. [5–7]). In contrast to the periodically driven system, investigation of random noise driven oscillators of system (1) has been paid little attention to. In Ref. [8], the present authors studied the effect of bounded noise on the chaotic motion of a harmonically excited elastic beam having three stable and two unstable equilibrium positions. By the random Melnikov method together with its associated mean-square criterion and the numerically calculated top Lyapunov exponents, the results showed that for larger noise intensity the threshold amplitude of bounded noise moving upwards with the increase of noise intensity, while contradiction appeared for smaller noise intensity. For the case of Eq. (1) with cubic nonlinearity only, Wei and Leng [9] have computed the top Lyapunov exponents of DuffingÕs equation perturbed by Gaussian white noise, leading to the conclusion that noise tends to stabilize the system. By Melnikov process, Lin and Yim [10] have studied the periodically forced Duffing system with external random perturbation, revealing that the presence of noise lowers the threshold and enlarges the possible chaotic domain in parameter space. Xie [11] discussed the effect of noise on the chaotic behavior of a buckled column under parametric excitation by a modified Melnikov method, and found that noise tending to increase the homoclinic threshold; meanwhile, Lyapunov exponents were simulated and the results revealed that the critical value of periodic forcing amplitude for onset of chaotic motion being decreased with the increase of noise intensity. Liu et al. [12] studied the effect of bounded noise on Duffing system with a homoclinic orbit under parametric excitation, and found that for larger noise intensity the threshold amplitude of bounded noise moving upward with the increase of noise intensity, while contradiction was found for smaller noise intensity. Therefore, the effect of random noise on chaotic behavior of Duffing oscillator is still riddling and further investigation is needed. The main purpose of this paper is to investigate the chaotic motion of the Duffing oscillator in a /6 potential under harmonic and Gaussian white noise excitations. The organization of this paper is as follows. In Section 2, theoretical analysis is presented, that is the derivation of necessary conditions for chaos as well as sufficient conditions for appearance of fractal basin boundaries due to homoclinic or heteroclinic bifurcation by virtue of Melnikov theory. In Section 3, numerical simulations are carried out to obtain the conditions of homoclinic and heteroclinic bifurcations, fractal basin boundaries, top Lyapunov exponents, Poincare maps of the system (1). Finally, we end with some concluding remarks.
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
129
2. Random Melnikov method Throughout the paper, we use the following values for the system parameters b = 1.0, c = 0.5 and d = 0.038. The unperturbed system, i.e. system (1) with a = c = r = 0, x_ ¼ y; ð2Þ y_ ¼ bx cx3 dx5 is a Hamiltonian system with Hamiltonian function 1 1 1 1 H ðx; yÞ ¼ y 2 þ bx2 þ cx4 þ dx6 . 2 2 4 6
ð3Þ
Through the analysis of the fixed points (xi, yi) and their stability for Eq. (2), one can see that there exits five fixed points: si(i = 1, 2) being saddles, Ci(i = 0, 1, 2) being centers as shown in Fig. 1(a). There is a homoclinic orbit C ho to þ and from S2 and there is another homoclinic orbit Cþ ho to and from S1. And there is a heteroclinic orbit Che from S1 to S2, together with another heteroclinic orbit Cþ he from S2 to S1. The phase portrait and the potential function of system (2) are shown in Fig. 1(a) and (b), respectively. Melnikov technique [4], proposed by Melnikov himself, is one of the few analytical methods available for determining the existence of the chaotic motion in near-integrable system subject to dissipative time-dependent perturbation. The main idea of Melnikov technique is to find a function that can measure the distance between the stable and unstable manifolds for a saddle or two saddles of the perturbed system. If the function vanishes for a certain bifurcation parameter value, then the stable and unstable manifolds will intersect each other away from the saddle point or points in the Poincare section. By a theorem attributed to Poincare [4], if the stable and unstable manifolds cross each other once, they will intersect an infinite number of times, thus forming a type of SmaleÕs horseshoes mapping leading to chaos. The Melnikov technique was firstly applied by Holmes [13] to study the chaotic attractor of a periodically driven Duffing oscillator with negative linear stiffness. A generalized version of the Melnikov function for a system subject to an excitation with multiple frequencies was introduced by Wiggins [14,15], and then Frey and Simiu [16] presented a generalized random Melnikov technique to study the effect of noise on near-integrable second-order dynamical system. Now supposed that a, c and r are small parameters with the same order as e, that is a ¼ e a, c ¼ ec, r ¼ e r, and (x0, y0) = (x0(t), y0(t)) is a homoclinic or heteroclinic orbit of the unperturbed system (2). Then the random Melnikov process can be obtained by using the formula given by Wiggins [14] Z þ1 Z þ1 Z þ1 cy 0 cos Xðt1 tÞ dt þ ay 20 dt þ y 0 nt2 t dt r Mðt1 ; t2 Þ ¼ ð4Þ 1 1 1 ¼ I þ zðt1 Þ þ Z t2 ; where the first two integrals in Eq. (4) represent the mean of the random Melnikov process due to damping and har y 0 ðtÞ monic excitation, and the last integral denotes the random portion due to random noise. The expression hðtÞ ¼ r can be regarded as the impulse response function of an invariant linear system while n(t) is an input of the system. Thus the variance of r2z as the output of the system can be obtained Z þ1 jH ðxÞj2 S n ðxÞ dx; ð5Þ r2z ¼ 1
where Sn(x) is the spectral density n(t) and H(x) is the frequency response function of the system obtained through R þ1 ofjxt y 0 e Fourier transform H ðxÞ ¼ 1 r dt. Since it is difficult to give analytical solution of y0(t) from Eq. (2), numerical calculation is employed to solve y0(t) in Section 3.1. Note that y0(t) is a function of time t from 1 to +1. It is seen from the phase portrait of system (2) displayed in Fig. 1(a) that it is convenient to choose a starting point p1 which is an intersecting point of the homoclinic þ orbit Cþ ho with x-axis, and a starting point p2 which is an intersecting point of heteroclinic orbit Che with y-axis, then y0(t) will be an odd function of time for the homoclinic orbit and even function for the heteroclinic orbit. Therefore, for the homoclinic orbit the Melnikov function can be simplified as M hom ðt1 ; t2 Þ ¼ 2aA þ 2c sin Xt1 I hom ðXÞ þ Z t2 ; ð6Þ R þ1 2 R þ1 where A ¼ 0 y 0 dt is a constant once y0 is given, and I hom ðXÞ ¼ 0 y 0 sin Xt dt is function of the frequency X. The criterion for possible chaotic motion based on Melnikov process is performed in mean-square representation
130
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
h2aAi2 ¼ h2c sin Xt1 I hom ðXÞi2 þ
ð2 rÞ2 p
Z 0
þ1
Z
0
þ1
2 y 0 sin xt dt dx.
ð7Þ
Then the mean-square criterion for possible chaotic motion in terms of parameters (c, X, a) of system (1) is given by Z 2 Z r2 þ1 þ1 dx. y sin xt dt ð8Þ haAi2 6 hcI hom ðXÞi2 þ 0 p 0 0 For the heteroclinic orbit the Melnikov function can be simplified as ð9Þ M het ðt1 ; t2 Þ ¼ 2aB þ 2c cos Xt1 I her ðXÞ þ Z t2 ; R þ1 2 R þ1 where B ¼ 0 y 0 dt is a constant once y0 is given, I her ðXÞ ¼ 0 y 0 cos Xt dt is function of the frequency X. The criterion for possible chaotic response based on Melnikov process is performed in mean-square representation Z 2 Z ð2 rÞ2 þ1 þ1 dx. h2aBi2 ¼ h2c cos Xt1 I her ðXÞi2 þ y cos xt dt ð10Þ 0 p 0 0 Then the mean-square criterion for possible chaotic response in terms of parameters (c, X, a) of system (1) is given by 2 Z Z r2 þ1 þ1 2 2 y 0 cos xt dt dx. ð11Þ haBi 6 hcI her ðXÞi þ p 0 0 When system (1) is free of noise, i.e. r = O in Eq. (2), system (1) degenerates into a Duffing oscillator under external
periodic excitation. For the homoclinic orbits C
ho and heteroclinic orbits Che , the Melnikov function can be simplified as M hom ðt1 Þ ¼ 2aA þ 2c sin Xt1 I hom ðXÞ;
ð12Þ
Fig. 2. Upper bound for possible chaotic domain due to homoclinic bifurcation: r = 0.0 and r = 1.0 for solid and dashed lines, respectively. (a) Upper surface in (c, X, a) space with r = 1.0; (b) upper surface in (c, X, a) space with r = 0.0; (c) upper bound in (c, a) plane with X = 1.0, and (d) upper bound in (X, a) plane with c = 0.8.
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
131
and M het ðt1 Þ ¼ 2aB þ 2c cos Xt1 I her ðXÞ;
ð13Þ
with respectively. Then in terms of parameters (c, X, a) of system (1), the possible chaotic response due to homoclinic and heterclinic bifurcation can be given by jaAj 6 jcI hom ðXÞj;
ð14Þ
jaBj 6 jcI her ðXÞj;
ð15Þ
and
respectively.
3. Numerical calculations 3.1. Simulation of the bifurcation conditions for chaos The algorithms, for calculating the bifurcation conditions due to homoclinic bifurcation or heterclinic bifurcation, essentially include two steps. Firstly, Eq. (2) is numerically solved for y0(t) using the fourth order Runge–Kutter technique with the time step 0.001, and then we numerically integrate the integrals A, B, Ihom and Iher with the time step 0.1. Even for those integrals infinite upper limits, only finite limits are needed in practical computation, since y0(t) exponentially converges to the corresponding saddle point at which y0(1) = 0. For system (2), the total length of time in
Fig. 3. Upper bound for possible chaotic domain due to heteroclinic bifurcation: r = 1.0 and r = 0.0 for solid and dashed lines, respectively. (a) Upper surface in (c, X, a) space with r = 1.0, (b) upper surface in (c, X, a) space with r = 0.0, (c) upper bound in (c, a) plane with X = 1.0, and (d) upper bound in (X, a) plane with c = 0.8.
132
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
computation is 10 non-dimensional time units for both the homoclinic and heteroclinic orbits. Once y0(t) is confirmed, the integrals A, B, Ihom and Iher are simply a matter of finite summations over the time interval 0–10.
Fig. 4. The threshold amplitude c versus noise intensity r for (a) homoclinic bifurcation, (b) heteroclinic bifurcation (—: analytic result, - -: numerical result).
Fig. 5. Basins of attraction for motion around x = x1. (a) c = 0.35, (b) c = 0.6258, and (c) c = 0.7.
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
133
The upper bounds of possible chaotic domain due to the homoclinic bifurcation for noise-free and noise system (1) are obtained by equating the expressions in Eqs. (8) and (14) and portrayed in Fig. 2(a) and (b), respectively. They depict the domains of possible occurrence of chaotic motion in parameters (c, X, a) space. The upper bounds are also presented by the dashed and solid lines for the cases with and without noise perturbation in (c, a) plane and in (X, a) plane as shown in Fig. 2(c) and (d), respectively. Similarly, the upper bounds of possible chaotic domain due to the heteroclinic bifurcation for noise-free and noise system (1) are obtained by equating the expressions in Eqs. (11) and (15) and portrayed in Fig. 3(a) and (b), respectively. They depict the domains of possible occurrence of chaotic motion in (c, X, a) space as that of homocilnic bifurcation stated above. The upper bounds are also presented by the dashed and solid lines for the cases with and without noise perturbation in (c, a) plane and in (X, a) plane as shown in Fig. 3(c) and (d), respectively. From Figs. 2 and 3 one can see that the presence of noise lowers the threshold and enlarges the possible chaotic domain in parameter space. The mean-square criteria in Eqs. (8) and (11) also define the threshold amplitude of harmonic excitation for onset of possible chaos in system (1), respectively. When the parameters are selected as a = 0.678, X = 1.0, the thresholds c versus noise intensity r are shown in Fig. 4.
Fig. 6. Basins of attraction for motion around x = x0. (a) c = 0.35, (b) c = 0.6258, and (c) c = 0.7.
134
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
3.2. Fractal basin boundaries Basin boundary simulations [6,7,17,18] are one means for checking homoclinic or heteroclinic bifurcation. The result of Li and Moon [6] implies that a homoclinic or heteroclinic bifurcation is a sufficient condition for the appearance of fractal basin boundaries. Such fractal boundaries indicate that whether the system is attracted to one or the other periodic attractor may be very sensitive to initial conditions. Thus small uncertainties in the initial conditions can lead to unpredictability of the system output even when the motion is not chaotic. Numerical simulations are carried out in order to determine the basins boundaries for noise-free system (1) for different values of the harmonic excitation amplitude c at a = 0.3, X = 1.0. The system has three attractors in the neighborhood of x = x0 and x = ± x1. For small c, the steady state motion will be periodic in one of the three attractors. Even if c is increased beyond the critical value for homoclinic bifurcation, it is still possible that the final steady motion could be periodic rather than chaotic. By performing a scan of the initial conditions in the (x, y) plane for various values of c, we find that when c is less than the homoclinic critical value 0.6255, the basins of attraction(marked regions) are regular(see Figs. 5(a), 6(a), and 7(a)). As c increases, the regular shape of basin of attraction is destroyed and the fractal behavior becomes more and more visible (see Figs. 5(b), (c), 6(b), (c), and 7(b), (c)). Note that the results of Fig. 5 represent the basins of attraction of the motion around the well x = x1. Figs. 6 and 7 correspond to the basins of attraction of the motion with different harmonic excitation values c around the well x = x0 and x = x1, respectively.
Fig. 7. Basins of attraction for motion around x = x1. (a) c = 0.35, (b) c = 0.6258, and (c) c = 0.7.
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
135
Fig. 8. The top Lyapunov exponent versus the harmonic excitation amplitude for noise-free system with (a) r = 0 and noisy system with (b) r = 0.01, (c) r = 0.1, (d) r = 0.3, (e) r = 0.5 and (f) r = 0.8.
136
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
3.3. The top Lyapunov exponent The Lyapunov exponents characterize quantitatively the dynamics of a system representing the asymptotic rate of exponential convergence or divergence of nearby orbits in phase space. Exponential divergence of nearby orbits implies that the dynamical behavior is sensitive to initial conditions. A positive top Lyapunov exponent for a bounded attractor is usually a sign of chaos. To check the threshold of harmonic excitation amplitude for onset of possible chaos obtained in Section 3.1, the top Lyapunov exponents of system (1) are also calculated by the Wolf et al. [19] algorithm. Numerical calculations have been made for the following parameter values: a = 0.678, X = 1.0 and the results for top Lyapunov exponents versus harmonic excitation amplitude of systems (1) are displayed in Fig. 8 for some different noise intensity values. Meanwhile, the thresholds of harmonic excitation amplitude for onset of chaos are obtained by letting the top Lyapunov exponents vanish, which are also shown in Fig. 4. From Fig. 8 one can see that for smaller values of harmonic excitation amplitude c, the top Lyapunov exponent is negative. As c increases, the top Lyapunov exponent changes from negative value to positive value, signifying the presence of chaotic motion. In the absence of noise, i.e. r = 0, beyond the threshold c for onset of chaotic motion, there are some ‘‘periodic windows’’, in which k becomes negative again and the system is then periodic. However, for larger noise intensity (r > rc, rc 0.098), there are no longer ‘‘periodic windows’’ presented up to the cutoff value of c in simulation. It is seen from Fig. 4 that the result obtained by Melnikov process with its associated mean-square criteria and that obtained by the numerical calculation of the top Lyapunov exponent yield the same variation trend, i.e. the threshold c decreases as noise intensity r increases, which implies that the larger the noise intensity results in the more possible chaotic domain in parameter space. 3.4. Poincare map To obtain a vivid picture about the behavior of system (1), we also investigate system (1) through its Poincare maps. The successive iteration of Poincare map is defined as,
Fig. 9. Poincare maps for noise-free with (a) r = 0.06, c = 1.6, and noisy system with (b) r = 0.08, c = 1.6, (c) r = 0.1, c = 1.6.
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
137
Fig. 10. Poincare maps for noise-free with (a) r = 0.0, c = 2.0, and noisy system with (b) r = 0.05, c = 2.0, (c) r = 0.1, c = 2.0.
P:
X
!
X ;
X
¼ fðxðtÞ; x_ ðtÞÞjt ¼ 0; 2p=X; 4p=X; . . .g R2 .
One hundred initial points are randomly chosen on the phase plane. For each initial condition, the differential Eq. (1) is solved by four order Rutter–Kutter method and the solution is plotted for every T = 2p/X. For each initial point, after deleting the first 1000 transients, the succeeded 19000 iteration points are plotted. It is known from analysis of the top Lyapunov exponents that the noise-free system (1) is periodic for c = 1.6 and chaotic for c = 2.0. Results are displayed in Fig. 9 for c = 1.6, r = 0.06, 0.08, 0.1 and Fig. 10 for c = 2.0, r = 0.0, 0.05, 0.1. It is seen from Fig. 9 that for larger noise intensity the map occupies larger area in phase plane. From Fig. 10 one can seen that the chaotic attractor is diffused by random noise, and the larger noise intensity results in the more diffused chaotic attractor.
4. Conclusion In this paper, the effect of Gaussian white noise on the chaotic motion of a harmonically excited Duffing oscillator in a triple well /6 potential possessing both homoclinic and heteroclinic orbits is detailed investigated. Based on homoclinic or heteroclinic bifurcation, necessary conditions for possible chaotic motion as well as sufficient conditions for the fractal basin boundaries are derived by Melnikov theory, which are employed to determine the fractality of basin of attraction, to depict the bifurcation surfaces in different parameter space, and to establish the threshold of harmonic excitation amplitude for onset of chaos. The results indicate that the presence of noise lowers the threshold amplitude and enlarges the possible chaotic domain in parameter space; moreover, the threshold of harmonic excitation amplitude for onset of chaos decreases as the intensity of noise increases. Numerical calculation of the top Lyapunov exponents of the original system also validates that the threshold amplitude for onset of chaos will move downwards as the noise intensity increases. Then one may conclude that the larger the noise intensity results in the more possible chaotic
138
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
domain in parameter space. The effect of noise on the system response is also investigated through its Poincare maps and the results indicate that the map occupies larger area in phase plane for larger noise intensity.
Acknowledgments This work was supported by the National Natural Science Foundation of China (Grant No. 10472091) and (Grant No. 10332030) and NSF of ShaanÕxi Province (Grant No. 2003A03).
References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19]
Thompson JM, Stewart HB. Non-linear dynamics and chaos. New York: Wiley; 1986. Ueda YJ. Stat Phys 1979;20:191. Moon FC. Chaotic and fractal dynamics. New York: Wiley; 1992. Guckenheimer J, Holmes P. Nonlinear oscillations, dynamical system, and bifurcations of vector fields. New York: SpringerVerlag; 1983. Zhujun J, Ruiqi W. Chaos, Solitions & Fractals 2005;22(3):887. Li GX, Moon FC. J Sound Vib 1990;136(1):17. Tchoukuegno R, Nana Bbendjio BR, Woafo P. Physica A 2002;34:362. Xiaoli Y, Wei X, Zhongkui S, Tang F. Chaos, Solitions & Fractals 2005;25:415. Wei JG, Leng G. Appl Math Comput 1997;88:77. Lin H, Yim SCS. J Appl Mech ASME 1996;63:509. Xie WC. Nonlinear Stochastic Dyn 1994;AMD-192(DE-73):215. Liu WY, Zhu WQ, Huang ZL. Chaos, Solitons & Fractals 2001;12:527. Holmes PJ. Philos Trans R Soc A 1979;292:419. Wiggins S. Global bifurcations and chaos: analytical methods. New York: Springer-Verlag; 1988. Wiggins S. Introduction to applied nonlinear dynamical system and chaos. New York: Springer-Verlag; 1990. Frey M, Simiu E. Physica D 1993;63:321. Moon FC, Li GX. Phys Rev Lett 1985;55:1439. Grebogi C, Macdonald S, Ott E, Yorke J. Phys Lett A 1983;99:415. Wolf A, Swift J, Swinney H, Vastano A. Physica D 1985;16:285.
Effect of random noise on chaotic motion of a particle in a /6 potential Zhongkui Sun a, Wei Xu a
a,*
, Xiaoli Yang
a,b
Department of Applied Mathematics, Northwestern Polytechnic University, Xi’an 710072, PR China b Department of Mathematics, Shaan’xi Normal University, Xi’an 710062, PR China Accepted 1 February 2005
Abstract The chaotic behaviors of a particle in a triple well /6 potential possessing both homoclinic and heteroclinic orbits under harmonic and Gaussian white noise excitations are discussed in detail. Following Melnikov theory, conditions for the existence of transverse intersection on the surface of homoclinic or heteroclinic orbits for triple potential well case are derived, which are complemented by the numerical simulations from which we show the bifurcation surfaces and the fractality of the basins of attraction. The results reveal that the threshold amplitude of harmonic excitation for onset of chaos will move downwards as the noise intensity increases, which is further verified by the top Lyapunov exponents of the original system. Thus the larger the noise intensity results in the more possible chaotic domain in parameter space. The effect of noise on Poincare maps is also investigated. Ó 2005 Published by Elsevier Ltd.
1. Introduction Collective investigations have been dedicated to the dynamics of a harmonically excited particle in a /4 potential [1–4] 1 1 V ðxÞ ¼ bx2 þ cx4 ; 2 4 where b and c are constants. Nevertheless, with the advent of the study of chaotic motion by means of strange attractors, top Lyapunov exponents, Poincare maps and fractal basin boundaries, it has become necessary to seek for better understanding of nonlinear system with higher order nonlinear terms. In this paper, we consider the extended Duffing oscillator by introducing a quintic term described by the following equation: x_ ¼ y; ð1Þ y_ ¼ bx cx3 dx5 ay þ c cos Xt þ rnðtÞ;
*
Corresponding author. Tel.: +86 29 88495453; fax: +86 29 88494174. E-mail addresses: [email protected] (Z. Sun), [email protected] (W. Xu).
0960-0779/$ - see front matter Ó 2005 Published by Elsevier Ltd. doi:10.1016/j.chaos.2005.02.033
128
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
Fig. 1. Phase portrait and potential of system (2) with b = 1.0, c = 0.5, d = 0.038. (a) Phase portrait and (b) potential.
where b, c and d are constants parameters; a is the damping ratio; n(t) represents a standard Gaussian white noise with Dirac-Delta function as its correction, namely, En(t) = 0, En(t)n(t + s) = d(s) and c and r denote the amplitude or level of the harmonic or noise excitation, respectively. The interest in such a mode is due to the fact that it is a universal nonlinear differential equation, and many nonlinear oscillators in physical, engineering and biological problems can be really described by the model or analogous ones. The potential from which is a /6 potential given by 1 1 1 V ðxÞ ¼ bx2 þ cx4 þ dx6 . 2 4 6 Depending on the set of the parameters, it can be considered at least three physically interesting situations where the potential is (i) single well, (ii) double well, and (iii) triple well. Throughout this paper, our analysis is carried out on the case of triple well with a double hump which does not lead to unbounded motion (see Fig. 1(b)). The effects of excitation, consisting of harmonic and random noise excitations, on nonlinear dynamical systems exhibiting chaotic behavior have been examined in recent years [5–12]. Especially the chaotic behaviors of this type of nonlinear system under periodic excitation have been extensively investigated (for instance, see Refs. [5–7]). In contrast to the periodically driven system, investigation of random noise driven oscillators of system (1) has been paid little attention to. In Ref. [8], the present authors studied the effect of bounded noise on the chaotic motion of a harmonically excited elastic beam having three stable and two unstable equilibrium positions. By the random Melnikov method together with its associated mean-square criterion and the numerically calculated top Lyapunov exponents, the results showed that for larger noise intensity the threshold amplitude of bounded noise moving upwards with the increase of noise intensity, while contradiction appeared for smaller noise intensity. For the case of Eq. (1) with cubic nonlinearity only, Wei and Leng [9] have computed the top Lyapunov exponents of DuffingÕs equation perturbed by Gaussian white noise, leading to the conclusion that noise tends to stabilize the system. By Melnikov process, Lin and Yim [10] have studied the periodically forced Duffing system with external random perturbation, revealing that the presence of noise lowers the threshold and enlarges the possible chaotic domain in parameter space. Xie [11] discussed the effect of noise on the chaotic behavior of a buckled column under parametric excitation by a modified Melnikov method, and found that noise tending to increase the homoclinic threshold; meanwhile, Lyapunov exponents were simulated and the results revealed that the critical value of periodic forcing amplitude for onset of chaotic motion being decreased with the increase of noise intensity. Liu et al. [12] studied the effect of bounded noise on Duffing system with a homoclinic orbit under parametric excitation, and found that for larger noise intensity the threshold amplitude of bounded noise moving upward with the increase of noise intensity, while contradiction was found for smaller noise intensity. Therefore, the effect of random noise on chaotic behavior of Duffing oscillator is still riddling and further investigation is needed. The main purpose of this paper is to investigate the chaotic motion of the Duffing oscillator in a /6 potential under harmonic and Gaussian white noise excitations. The organization of this paper is as follows. In Section 2, theoretical analysis is presented, that is the derivation of necessary conditions for chaos as well as sufficient conditions for appearance of fractal basin boundaries due to homoclinic or heteroclinic bifurcation by virtue of Melnikov theory. In Section 3, numerical simulations are carried out to obtain the conditions of homoclinic and heteroclinic bifurcations, fractal basin boundaries, top Lyapunov exponents, Poincare maps of the system (1). Finally, we end with some concluding remarks.
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
129
2. Random Melnikov method Throughout the paper, we use the following values for the system parameters b = 1.0, c = 0.5 and d = 0.038. The unperturbed system, i.e. system (1) with a = c = r = 0, x_ ¼ y; ð2Þ y_ ¼ bx cx3 dx5 is a Hamiltonian system with Hamiltonian function 1 1 1 1 H ðx; yÞ ¼ y 2 þ bx2 þ cx4 þ dx6 . 2 2 4 6
ð3Þ
Through the analysis of the fixed points (xi, yi) and their stability for Eq. (2), one can see that there exits five fixed points: si(i = 1, 2) being saddles, Ci(i = 0, 1, 2) being centers as shown in Fig. 1(a). There is a homoclinic orbit C ho to þ and from S2 and there is another homoclinic orbit Cþ ho to and from S1. And there is a heteroclinic orbit Che from S1 to S2, together with another heteroclinic orbit Cþ he from S2 to S1. The phase portrait and the potential function of system (2) are shown in Fig. 1(a) and (b), respectively. Melnikov technique [4], proposed by Melnikov himself, is one of the few analytical methods available for determining the existence of the chaotic motion in near-integrable system subject to dissipative time-dependent perturbation. The main idea of Melnikov technique is to find a function that can measure the distance between the stable and unstable manifolds for a saddle or two saddles of the perturbed system. If the function vanishes for a certain bifurcation parameter value, then the stable and unstable manifolds will intersect each other away from the saddle point or points in the Poincare section. By a theorem attributed to Poincare [4], if the stable and unstable manifolds cross each other once, they will intersect an infinite number of times, thus forming a type of SmaleÕs horseshoes mapping leading to chaos. The Melnikov technique was firstly applied by Holmes [13] to study the chaotic attractor of a periodically driven Duffing oscillator with negative linear stiffness. A generalized version of the Melnikov function for a system subject to an excitation with multiple frequencies was introduced by Wiggins [14,15], and then Frey and Simiu [16] presented a generalized random Melnikov technique to study the effect of noise on near-integrable second-order dynamical system. Now supposed that a, c and r are small parameters with the same order as e, that is a ¼ e a, c ¼ ec, r ¼ e r, and (x0, y0) = (x0(t), y0(t)) is a homoclinic or heteroclinic orbit of the unperturbed system (2). Then the random Melnikov process can be obtained by using the formula given by Wiggins [14] Z þ1 Z þ1 Z þ1 cy 0 cos Xðt1 tÞ dt þ ay 20 dt þ y 0 nt2 t dt r Mðt1 ; t2 Þ ¼ ð4Þ 1 1 1 ¼ I þ zðt1 Þ þ Z t2 ; where the first two integrals in Eq. (4) represent the mean of the random Melnikov process due to damping and har y 0 ðtÞ monic excitation, and the last integral denotes the random portion due to random noise. The expression hðtÞ ¼ r can be regarded as the impulse response function of an invariant linear system while n(t) is an input of the system. Thus the variance of r2z as the output of the system can be obtained Z þ1 jH ðxÞj2 S n ðxÞ dx; ð5Þ r2z ¼ 1
where Sn(x) is the spectral density n(t) and H(x) is the frequency response function of the system obtained through R þ1 ofjxt y 0 e Fourier transform H ðxÞ ¼ 1 r dt. Since it is difficult to give analytical solution of y0(t) from Eq. (2), numerical calculation is employed to solve y0(t) in Section 3.1. Note that y0(t) is a function of time t from 1 to +1. It is seen from the phase portrait of system (2) displayed in Fig. 1(a) that it is convenient to choose a starting point p1 which is an intersecting point of the homoclinic þ orbit Cþ ho with x-axis, and a starting point p2 which is an intersecting point of heteroclinic orbit Che with y-axis, then y0(t) will be an odd function of time for the homoclinic orbit and even function for the heteroclinic orbit. Therefore, for the homoclinic orbit the Melnikov function can be simplified as M hom ðt1 ; t2 Þ ¼ 2aA þ 2c sin Xt1 I hom ðXÞ þ Z t2 ; ð6Þ R þ1 2 R þ1 where A ¼ 0 y 0 dt is a constant once y0 is given, and I hom ðXÞ ¼ 0 y 0 sin Xt dt is function of the frequency X. The criterion for possible chaotic motion based on Melnikov process is performed in mean-square representation
130
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
h2aAi2 ¼ h2c sin Xt1 I hom ðXÞi2 þ
ð2 rÞ2 p
Z 0
þ1
Z
0
þ1
2 y 0 sin xt dt dx.
ð7Þ
Then the mean-square criterion for possible chaotic motion in terms of parameters (c, X, a) of system (1) is given by Z 2 Z r2 þ1 þ1 dx. y sin xt dt ð8Þ haAi2 6 hcI hom ðXÞi2 þ 0 p 0 0 For the heteroclinic orbit the Melnikov function can be simplified as ð9Þ M het ðt1 ; t2 Þ ¼ 2aB þ 2c cos Xt1 I her ðXÞ þ Z t2 ; R þ1 2 R þ1 where B ¼ 0 y 0 dt is a constant once y0 is given, I her ðXÞ ¼ 0 y 0 cos Xt dt is function of the frequency X. The criterion for possible chaotic response based on Melnikov process is performed in mean-square representation Z 2 Z ð2 rÞ2 þ1 þ1 dx. h2aBi2 ¼ h2c cos Xt1 I her ðXÞi2 þ y cos xt dt ð10Þ 0 p 0 0 Then the mean-square criterion for possible chaotic response in terms of parameters (c, X, a) of system (1) is given by 2 Z Z r2 þ1 þ1 2 2 y 0 cos xt dt dx. ð11Þ haBi 6 hcI her ðXÞi þ p 0 0 When system (1) is free of noise, i.e. r = O in Eq. (2), system (1) degenerates into a Duffing oscillator under external
periodic excitation. For the homoclinic orbits C
ho and heteroclinic orbits Che , the Melnikov function can be simplified as M hom ðt1 Þ ¼ 2aA þ 2c sin Xt1 I hom ðXÞ;
ð12Þ
Fig. 2. Upper bound for possible chaotic domain due to homoclinic bifurcation: r = 0.0 and r = 1.0 for solid and dashed lines, respectively. (a) Upper surface in (c, X, a) space with r = 1.0; (b) upper surface in (c, X, a) space with r = 0.0; (c) upper bound in (c, a) plane with X = 1.0, and (d) upper bound in (X, a) plane with c = 0.8.
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
131
and M het ðt1 Þ ¼ 2aB þ 2c cos Xt1 I her ðXÞ;
ð13Þ
with respectively. Then in terms of parameters (c, X, a) of system (1), the possible chaotic response due to homoclinic and heterclinic bifurcation can be given by jaAj 6 jcI hom ðXÞj;
ð14Þ
jaBj 6 jcI her ðXÞj;
ð15Þ
and
respectively.
3. Numerical calculations 3.1. Simulation of the bifurcation conditions for chaos The algorithms, for calculating the bifurcation conditions due to homoclinic bifurcation or heterclinic bifurcation, essentially include two steps. Firstly, Eq. (2) is numerically solved for y0(t) using the fourth order Runge–Kutter technique with the time step 0.001, and then we numerically integrate the integrals A, B, Ihom and Iher with the time step 0.1. Even for those integrals infinite upper limits, only finite limits are needed in practical computation, since y0(t) exponentially converges to the corresponding saddle point at which y0(1) = 0. For system (2), the total length of time in
Fig. 3. Upper bound for possible chaotic domain due to heteroclinic bifurcation: r = 1.0 and r = 0.0 for solid and dashed lines, respectively. (a) Upper surface in (c, X, a) space with r = 1.0, (b) upper surface in (c, X, a) space with r = 0.0, (c) upper bound in (c, a) plane with X = 1.0, and (d) upper bound in (X, a) plane with c = 0.8.
132
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
computation is 10 non-dimensional time units for both the homoclinic and heteroclinic orbits. Once y0(t) is confirmed, the integrals A, B, Ihom and Iher are simply a matter of finite summations over the time interval 0–10.
Fig. 4. The threshold amplitude c versus noise intensity r for (a) homoclinic bifurcation, (b) heteroclinic bifurcation (—: analytic result, - -: numerical result).
Fig. 5. Basins of attraction for motion around x = x1. (a) c = 0.35, (b) c = 0.6258, and (c) c = 0.7.
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
133
The upper bounds of possible chaotic domain due to the homoclinic bifurcation for noise-free and noise system (1) are obtained by equating the expressions in Eqs. (8) and (14) and portrayed in Fig. 2(a) and (b), respectively. They depict the domains of possible occurrence of chaotic motion in parameters (c, X, a) space. The upper bounds are also presented by the dashed and solid lines for the cases with and without noise perturbation in (c, a) plane and in (X, a) plane as shown in Fig. 2(c) and (d), respectively. Similarly, the upper bounds of possible chaotic domain due to the heteroclinic bifurcation for noise-free and noise system (1) are obtained by equating the expressions in Eqs. (11) and (15) and portrayed in Fig. 3(a) and (b), respectively. They depict the domains of possible occurrence of chaotic motion in (c, X, a) space as that of homocilnic bifurcation stated above. The upper bounds are also presented by the dashed and solid lines for the cases with and without noise perturbation in (c, a) plane and in (X, a) plane as shown in Fig. 3(c) and (d), respectively. From Figs. 2 and 3 one can see that the presence of noise lowers the threshold and enlarges the possible chaotic domain in parameter space. The mean-square criteria in Eqs. (8) and (11) also define the threshold amplitude of harmonic excitation for onset of possible chaos in system (1), respectively. When the parameters are selected as a = 0.678, X = 1.0, the thresholds c versus noise intensity r are shown in Fig. 4.
Fig. 6. Basins of attraction for motion around x = x0. (a) c = 0.35, (b) c = 0.6258, and (c) c = 0.7.
134
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
3.2. Fractal basin boundaries Basin boundary simulations [6,7,17,18] are one means for checking homoclinic or heteroclinic bifurcation. The result of Li and Moon [6] implies that a homoclinic or heteroclinic bifurcation is a sufficient condition for the appearance of fractal basin boundaries. Such fractal boundaries indicate that whether the system is attracted to one or the other periodic attractor may be very sensitive to initial conditions. Thus small uncertainties in the initial conditions can lead to unpredictability of the system output even when the motion is not chaotic. Numerical simulations are carried out in order to determine the basins boundaries for noise-free system (1) for different values of the harmonic excitation amplitude c at a = 0.3, X = 1.0. The system has three attractors in the neighborhood of x = x0 and x = ± x1. For small c, the steady state motion will be periodic in one of the three attractors. Even if c is increased beyond the critical value for homoclinic bifurcation, it is still possible that the final steady motion could be periodic rather than chaotic. By performing a scan of the initial conditions in the (x, y) plane for various values of c, we find that when c is less than the homoclinic critical value 0.6255, the basins of attraction(marked regions) are regular(see Figs. 5(a), 6(a), and 7(a)). As c increases, the regular shape of basin of attraction is destroyed and the fractal behavior becomes more and more visible (see Figs. 5(b), (c), 6(b), (c), and 7(b), (c)). Note that the results of Fig. 5 represent the basins of attraction of the motion around the well x = x1. Figs. 6 and 7 correspond to the basins of attraction of the motion with different harmonic excitation values c around the well x = x0 and x = x1, respectively.
Fig. 7. Basins of attraction for motion around x = x1. (a) c = 0.35, (b) c = 0.6258, and (c) c = 0.7.
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
135
Fig. 8. The top Lyapunov exponent versus the harmonic excitation amplitude for noise-free system with (a) r = 0 and noisy system with (b) r = 0.01, (c) r = 0.1, (d) r = 0.3, (e) r = 0.5 and (f) r = 0.8.
136
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
3.3. The top Lyapunov exponent The Lyapunov exponents characterize quantitatively the dynamics of a system representing the asymptotic rate of exponential convergence or divergence of nearby orbits in phase space. Exponential divergence of nearby orbits implies that the dynamical behavior is sensitive to initial conditions. A positive top Lyapunov exponent for a bounded attractor is usually a sign of chaos. To check the threshold of harmonic excitation amplitude for onset of possible chaos obtained in Section 3.1, the top Lyapunov exponents of system (1) are also calculated by the Wolf et al. [19] algorithm. Numerical calculations have been made for the following parameter values: a = 0.678, X = 1.0 and the results for top Lyapunov exponents versus harmonic excitation amplitude of systems (1) are displayed in Fig. 8 for some different noise intensity values. Meanwhile, the thresholds of harmonic excitation amplitude for onset of chaos are obtained by letting the top Lyapunov exponents vanish, which are also shown in Fig. 4. From Fig. 8 one can see that for smaller values of harmonic excitation amplitude c, the top Lyapunov exponent is negative. As c increases, the top Lyapunov exponent changes from negative value to positive value, signifying the presence of chaotic motion. In the absence of noise, i.e. r = 0, beyond the threshold c for onset of chaotic motion, there are some ‘‘periodic windows’’, in which k becomes negative again and the system is then periodic. However, for larger noise intensity (r > rc, rc 0.098), there are no longer ‘‘periodic windows’’ presented up to the cutoff value of c in simulation. It is seen from Fig. 4 that the result obtained by Melnikov process with its associated mean-square criteria and that obtained by the numerical calculation of the top Lyapunov exponent yield the same variation trend, i.e. the threshold c decreases as noise intensity r increases, which implies that the larger the noise intensity results in the more possible chaotic domain in parameter space. 3.4. Poincare map To obtain a vivid picture about the behavior of system (1), we also investigate system (1) through its Poincare maps. The successive iteration of Poincare map is defined as,
Fig. 9. Poincare maps for noise-free with (a) r = 0.06, c = 1.6, and noisy system with (b) r = 0.08, c = 1.6, (c) r = 0.1, c = 1.6.
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
137
Fig. 10. Poincare maps for noise-free with (a) r = 0.0, c = 2.0, and noisy system with (b) r = 0.05, c = 2.0, (c) r = 0.1, c = 2.0.
P:
X
!
X ;
X
¼ fðxðtÞ; x_ ðtÞÞjt ¼ 0; 2p=X; 4p=X; . . .g R2 .
One hundred initial points are randomly chosen on the phase plane. For each initial condition, the differential Eq. (1) is solved by four order Rutter–Kutter method and the solution is plotted for every T = 2p/X. For each initial point, after deleting the first 1000 transients, the succeeded 19000 iteration points are plotted. It is known from analysis of the top Lyapunov exponents that the noise-free system (1) is periodic for c = 1.6 and chaotic for c = 2.0. Results are displayed in Fig. 9 for c = 1.6, r = 0.06, 0.08, 0.1 and Fig. 10 for c = 2.0, r = 0.0, 0.05, 0.1. It is seen from Fig. 9 that for larger noise intensity the map occupies larger area in phase plane. From Fig. 10 one can seen that the chaotic attractor is diffused by random noise, and the larger noise intensity results in the more diffused chaotic attractor.
4. Conclusion In this paper, the effect of Gaussian white noise on the chaotic motion of a harmonically excited Duffing oscillator in a triple well /6 potential possessing both homoclinic and heteroclinic orbits is detailed investigated. Based on homoclinic or heteroclinic bifurcation, necessary conditions for possible chaotic motion as well as sufficient conditions for the fractal basin boundaries are derived by Melnikov theory, which are employed to determine the fractality of basin of attraction, to depict the bifurcation surfaces in different parameter space, and to establish the threshold of harmonic excitation amplitude for onset of chaos. The results indicate that the presence of noise lowers the threshold amplitude and enlarges the possible chaotic domain in parameter space; moreover, the threshold of harmonic excitation amplitude for onset of chaos decreases as the intensity of noise increases. Numerical calculation of the top Lyapunov exponents of the original system also validates that the threshold amplitude for onset of chaos will move downwards as the noise intensity increases. Then one may conclude that the larger the noise intensity results in the more possible chaotic
138
Z. Sun et al. / Chaos, Solitons and Fractals 27 (2006) 127–138
domain in parameter space. The effect of noise on the system response is also investigated through its Poincare maps and the results indicate that the map occupies larger area in phase plane for larger noise intensity.
Acknowledgments This work was supported by the National Natural Science Foundation of China (Grant No. 10472091) and (Grant No. 10332030) and NSF of ShaanÕxi Province (Grant No. 2003A03).
References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19]
Thompson JM, Stewart HB. Non-linear dynamics and chaos. New York: Wiley; 1986. Ueda YJ. Stat Phys 1979;20:191. Moon FC. Chaotic and fractal dynamics. New York: Wiley; 1992. Guckenheimer J, Holmes P. Nonlinear oscillations, dynamical system, and bifurcations of vector fields. New York: SpringerVerlag; 1983. Zhujun J, Ruiqi W. Chaos, Solitions & Fractals 2005;22(3):887. Li GX, Moon FC. J Sound Vib 1990;136(1):17. Tchoukuegno R, Nana Bbendjio BR, Woafo P. Physica A 2002;34:362. Xiaoli Y, Wei X, Zhongkui S, Tang F. Chaos, Solitions & Fractals 2005;25:415. Wei JG, Leng G. Appl Math Comput 1997;88:77. Lin H, Yim SCS. J Appl Mech ASME 1996;63:509. Xie WC. Nonlinear Stochastic Dyn 1994;AMD-192(DE-73):215. Liu WY, Zhu WQ, Huang ZL. Chaos, Solitons & Fractals 2001;12:527. Holmes PJ. Philos Trans R Soc A 1979;292:419. Wiggins S. Global bifurcations and chaos: analytical methods. New York: Springer-Verlag; 1988. Wiggins S. Introduction to applied nonlinear dynamical system and chaos. New York: Springer-Verlag; 1990. Frey M, Simiu E. Physica D 1993;63:321. Moon FC, Li GX. Phys Rev Lett 1985;55:1439. Grebogi C, Macdonald S, Ott E, Yorke J. Phys Lett A 1983;99:415. Wolf A, Swift J, Swinney H, Vastano A. Physica D 1985;16:285.