- Email: [email protected]
Crises and transient chaos in a driven Josephson junction
PhysicsLettersA 179 (1993) 337—342 North-Holland
PHYSICS LETTERS A
Crises and transient chaos in a driven Josephson junction F. Palmero Departamento de Fisica Aplicada, Universidad de Sevilla, P.O. Box 1065, Seville, Spain
and F. Romero Romero Departamento de Fjsjca Teorjca, Facultad deFIsica, Universidad deSevilla, P.O. Box 1065, Seville, Spain Received 26 January 1993; revised manuscript received 10 June 1993; accepted forpublication 14 June 1993 Communicatedby A.P. Fordy
Nonlinear dynamics of a single Josephson junction under radio-frequency excitation is analyzed by a model which takes into account interference effects and by means of numerical simulations. The study emphasizes the existence of a crisis, which is followedby chaotic transients over a portion of its state space. The mean lifetime is found to scale as I — I ~, where w is the frequency ofthe external excitation, a~,is its value at the crisis, and ö is the critical exponent.
1. Introduction Interest in the study of Josephson junctions has increased after the observation of Huberman et al. [1] that chaotic solutions might be responsible for the large noise increase in the Josephson parametric amplifiers [2]. Its mechanical analogue, the damped driven pendulum, has been the subject of extensive studies [3—5]and now it is clear that the Josephson junction is one of the powerful paradigms in nonlinear physics and that Josephson junction devices are imbued with a very rich dynamics and a variety of strange phenomena [6—10]. A typical Josephson junction can be described by McCumber’s model [11] and with a change of variables it is equivalent to the classical equation of motion of a simple pendulum with viscous damping, driven by a periodic external force. Nevertheless, to better describe the dynamics of the junction under the excitation by a radio-frequency current, the equation of motion must include a term that represents the interference effects between the quasiparticle-pair and quasiparticle currents, which is usually called the “cos 0” term. Most works about driven Josephson tunnelingjunctions have neglected
this interference term. The existing theories predict relative signs of the quasiparticle and quasiparticlepair currents which are in contradiction with some experimental results [12]. Then it appears reasonable to include in the evolution equations a “cos 0” term with a sign supported by experimental results. The analysis ofthe resulting equations might give rise to the prediction ofqualitatively new phenomena and might quantitatively affect the dynamical and fluctuation behaviour of Josephson devices. In this paper we have obtained computer solutions of the equation of motion modified by interference and analyzed the resulting data through various diagnostics tools: Poincaré sections, phase portraits, power spectra analysis, calculations of Lyapunov exponents and determination of basins of attraction. As the only varying control parameter we have considered the angular frequency of the external excitation. Thus, we have determined global bifurcation diagrams both for the noninterference model and for the interference model, which show substantial differences between them. By using the interpolated cell-mapping method (1CM) [131, we have obtained some basins of attraction, helping these to identify some topological
0375-9601/93/S 06.00 © 1993 Elsevier Science Publishers B.V. All rights reserved.
337
Volume 179, number 4,5
PHYSICS LETTERS A
types ofbifurcation. In this paper we are particularly interested in an interval of the control parameter where the system exhibits a chaotic blue sky catastrophe, associated with a hysteresis phenomenon, followed by transient chaos.
2. Model equation and numerical method McCumber’s model for a current-fed Josephson junction taking into account interference effects [141 has the following equivalent circuit equation, ~ + 4IEeR h (1 + y cos 0) 4xe dt2 0 = ‘o
sin (COrft)
+ I~sin 0
(1)
.
where 0 denotes the total phase difference across the junction, Cthe junction capacitor, h the Planck constant, e the electron charge, 4 the critical supercurrent, R0 the junction resistor corresponding to J~=ü, and y is the coherence coefficient which to 1o, of is therelated external the frequency, Wrt, and amplitude, excitation [15—171. This evolution equation can be transformed into the dimensionless form d2Ø ~
+a(l+ycosØ)~ +sinO=iesin(wr),
(2)
with T=(47~eI~/hC)”2t,a=R~’(h/4EeI~C)”2,~e= ‘o/j~~~ and w = w~/wo;wo = (4xeI~/hC)’/2~ In this form our equation has four control parameters, a, y, i~ and w, and can equivalently be written as
16 August 1993
considered values in the range 0.45 ~ w~0.73 and for ~e the fixed value ~e= 0.8. Some experiences on the Josephson plasma resonance in Pb—Pb tunnel junctions by Pedersen et al. [15] have shown a negative value for y, i.e. y= —0.9 ±0.2, then for our calculations we have chosen the value y= —0.8, which is in the range supported by the experimental results. Equations (3), (4) have been resolved in double precision arithmetic using a fourth-order variablestep-sizeRunge—Kutta method. In all cases, in order to eliminate transients, our solutions were run through over 200 periods of the driven current before actual data taking was started. For the numerical computations of the bifurcation diagrams the control parameter w has been increased (or decreased) in small steps. For each value of co we have obtained the steady-state solutions corresponding to a set of initial conditions randomly distributed in a determined subset of the phase space (X1, X2), and in a sufficient number such that all possible attractors are represented. Many aspects of the discussion are concerned with basins of attraction. Generally the determination of basins by numerical integration becomes extremely costly and time consuming. This approach has been used extensively in the literature when determining the basins of attraction of nonlinear systems [19]. Nevertheless, using the 1CM method one can determine the continuous time trajectory by spatially discretizing the system employing an interpolation procedure, with the advantage that one can generate an infinite number of distinct orbits from an Mx M mesh. In this paper we have used a 100 x 100 mesh over the area of phase space containing the attractors
dX 1 dr dX2 dr —
(3) =
—
a (1 + y cos X1 )X2
—
sin X1 + ~e sin (on), (4)
where the variable X2 is proportional to the difference of potential between the two superconductors. In order to obtain numerical solutions of eqs. (3), (4) we must choose a set of values for the parameters. Typical values of the coefficient a lie in the range 0.03 ~ a ~ 0.7, indicating that these systems are underdamped in the linear regime [11,18], and we have taken the fixed value a = 0.4. For w we have 338
of interest. For our model this grid is adequate to generate all the attractors that can be found by numerical integration, and the corresponding basins of attraction have enough resolution.
3. Bifurcation diagrams and crisis In order to appreciate to some extent the dynamical behaviour of the interference model as the parameter Co is varied, we have determined the global bifurcation diagram corresponding to the variable X2. Also, we have determined the global bifurcation diagrams for the corresponding noninterference model,
Volume 179, number 4,5
PHYSICS LETTERS A
16 August 1993
i.e. with the same values for all the parameters except for y, which now is 0. Figures 1 and 2 show both diagrams. From the comparison it is evident that the interference term modifies the Poincaré map rather substantially. Many aspects concerning the dynamics of the interference model could be analyzed, but
By decreasing the values of w from 0.73 we find that at w~0.672 a saddle node bifurcation occurs originating a stable period-three orbit and simultaneously an unstable orbit. This period-three solution is symmetric, i.e. it has inversion symmetry about the point (—It, 0) in the plane (X1, X2). Equation
in this paper we want to emphasize only some phenomena occurring in an interval around w=0.666.
(2) is invariant under the transformation X—~’ X, r—’r+ (2n+ 1)x/w, where n is an integer, then if X(r) is a solution of eq. (2), so is —X(t+ (2n+ 1)x/co). For the periodic solution encountered these two solutions are essentially the same, i e they differ by an integer number ofcomplete revolutions, so that
________________________________________ ~
.~
~
___,
—
(5)
—X(r+(2n+l)x/@)+2n=X(t)
As the period-three orbit is created there exist also chaotic solutions, these were originated after a cascade of period doubling bifurcations, i.e. a Feigenbaum route to chaos occurred for larger values of cv.
-
4
>~
~
•
..~.
~
‘~
~•‘.
At co= 0.666 a catastrophic event is about to appear. With extensive use of the 1CM method we have determined the basins of attraction corresponding to these coexisting attractors If we supenmpose the corresponding bidimensional Poincare sections and that ofthe unstable orbit, we obtain a composite por-
‘~ \
45
5
55
6
/
~ )-~
-.
65
7
Fig. 1. Global bifurcation diagram for the variable X2 correspondmg to the interference model with a=0 4 y= —0 8 i =0 8 andwintherange0.45~o~0.73.
.
~
‘
•~:~:‘
—2
—1
~
0
1
2
4’
3
Fig. 3. Basins of attraction in the plane (X1, X2) for the chaotic
I
.5
‘
xi
:
.45
~ ‘
,,
—3
•~frI
.
~
.55
I
,.,,
.6
.65
.7
w Fig. 2. Global bifurcation diagram for the variable X2 corresponding to the noninterference model, i.e. y=0, with a=0.4, —0.8, i,=0.8 and coin the range 0.45 ~ w~0.71.
solutions (blank) and the coexistent period-three solution (shaded) computed from eqs. (3), (4) for w=0.666, a=0.4, y=—0.8andi~=0.8.A1soshownarethePoincar~sectionsofthe chaotic attractor, the period-three orbit (centers of the full cirdes) and the unstable orbit (center of the open circle). The patterned-dot structure forced by the pixel separation ofthe printer must be ignored, since it is not caused by the numerics.
339
Volume 179, number 4,5
PHYSICS LETTERS A
16 August 1993
trait which can be observed in fig. 3. The basin boundary passes through the open circles, i.e. it belongs to the stable manifold of the unstable orbit, and we can see that it is about to collide with the chaotic attractor. In fact, at co= w,.y~0.66575, a tangency of the stable and the unstable manifold ofa period-three orbit, located on the basin boundary and accessible from the basin boundary ofthe chaotic attractor, occurs, originating a discontinuous disappearance of the chaotic attractor or crisis [20]. Thompson and Stewart [21] reported the existence of a similar phenomenon found in the forced oscillator equation d2X —0.25 dX —X+X3=A sin(w’r). (6)
terference model (eq. (2)), and in the range of cv where the period-three orbit coexists with the chaotic attractor, the basin boundary has also a fractallike structure, as can be seen in fig. 3. After the catastrophe, the basin of attraction of the chaotic attractor is incorporated to the basin of attraction of the period-three orbit, and when the symmetrybreakingbifurcation occurs the basin of attraction is divided in two regions, each one corresponding to a different orbit, yielding the structure shown in fig. 4. In the interpretation ofthis figure we must take into account that the pixel separation of the printer detennines a structure not caused by the numerics. Thus, the large ligthly shaded gray area is not an extension only of the densely shaded basin, but also includes “blank fragments”, and the resulting basin
This crisis is closely related to a hysteresis phenomenon observed by Huberman and Crutchfield [22]. For our case the hysteresis phenomenon is a necessary condition. When decreasing cv, destruction of the chaotic attractor is made possible by a previous creation of a saddle orbit, but as this is accompanied by the creation of a stable orbit, there exist simultaneously two attractors, the chaotic attractor and the period-three orbit. A symmetry-breaking bifurcation of the stable period-three orbit occurs for cv very close to w~,originating two different asymmetric orbits, each one being an inversion of the other with respect to the point (—x, 0) in the plane
boundary is also fractal. These facts determine the high sensitivity to initial conditions, not only when there exists a chaotic attractor, but when the only two coexistent attractors are periodic. Recent works relative to driven oscillators with cubic and quartic potential wells [27,28] showed some comparable forms of basin-breakup transitions.
-~—
I
I
_______
(X 1, X2).
.
4. Basins of attraction Another important question to consider is the one about the structure of the boundaries between basins of attraction when there is coexistence of attractors. Thus, Yorke and collaborators in a series of papers [19,23—25] treat the problem of how basins of attraction change as a system parameter varies continuously, and which are the qualitative changes of the basin boundaries as a parameter passes through certam critical values. They called such changes “metamorphoses”, and presented numerical evidences for fractal boundaries between basins of attraction for nonchaotic attractors. Fractal boundaries appear in many nonlinear systems as for example in the noninterference McCumber model [26]. For the in340
.-
________
_________________ _____
_____
.
.
:
-
Fig. 4. Basins ofattraction in the plane (X1, X2) for the two coexistent penod-three orbit for w = 0.664. The large lightly shaded area corresponds to the interspersion of the blank area and the densely shaded area to the basins ofattraction ofthe two coexistent orbits.
Volume 179, number 4,5
PHYSICS LETTERS A
5. Transient chaos
16 August 1993
~
...,. ~
Once the chaotic attractor has disappeared and for cv~w~,every trajectory starting from a given initial condition ultimately comes to one of the two coex istent period three attractors Nevertheless, the times employed to reach the final states differ and strongly depend on the initial conditions specially if they pertain to the basin of attraction of the destroyed chaotic attractor i e the system exhibits transient chaos This type of phenomenom has been observed in many different scenanos, as for example in the parametnc damped pendulum [29], where chaotic transients arise from the spatial congruence of the basinofattractionofastablestationarysolutionand the remnant of the destroyed chaotic attractor. Grebogi et al. [30] also showed some phenomena of transient chaos associated to bidimensional maps that appear both after a boundary crisis and after an interior crisis. For these cases the average lifetime Tof a chaotic transient scales as I p I where p denotes the control parameter and p~its critical value, and ô is the critical exponent ofthe chaotic transient. In order to find the distribution of the chaotic transient times for a given value of cv an appropiate set of initial conditions must be selected. As a first step we have taken 40000 initial conditions randomly distributed in the rectangle of the phase space _It~Xi~It,—2.5i~X2~l.5.By using the 1CM —
‘i~.,
~
I ..
,
‘~
...,,.
~
‘~.
~ .
—. “~
~
3
~J
~
~, ~.,
~.
-2
I I
1
.
0
1
2
3
Xi
Fig. 5. Union ofthe basin cores of the two coexistent attractors for w=0.665. ‘II’
—
— ~,
method we have followed their trajectories and selected all the initial points which reach a determined final state after a time of less than five periods ofthe external driving force. This selection gives rise to what is called “the basin core of the attractor”. For cv~0.665 there are two coexistent periodic attractors, the union of the two basin cores for cv = 0.665 can be observed in fig. 5 and as it can be seen is very similar to the basin of attraction ofthe periodic orbit previous to the chaotic catastrophe. To study transient chaos the appropiate set of initial conditions is determined by eliminating from the set of 40000 points all of them which pertain to the basin core. For each selected point we have determined the number of cycles required for the system to evolve to the final states. For each value of cv we have determined a histogram illustrating the number of initial conditions N versus the number of cycles of transient duration NT.
—
~
s
—
—
20
-
__________________________________ 40
60
80
100
120
N1
Fig. 6. Histogram of grouped data obtained for w= 0.665.
In fig. 6 we show the histogram corresponding to cv = 0.66 5. Inspection of this figure suggests that the probability distribution P( t) for the chaotic transient lifetime T is of the form P( ‘r) = (CI T) exp ( T), then T represents the mean lifetime of the chaotic transient. A plot of ln (N) versus ln (Ni) shows a straight line with slope T~. In order to determine the type of dependence of Ton cv in the proximity of cv,~,we have calculated the mean lifetime T for some different values of cv. Representing log10 ( T) versus log10(cv— cv,,), as it is shown in fig. 7, we can —
—
341
Volume 179, number 4,5
PHYSICS LETTERS A
16 August 1993
References
~
[1]B.A. Huberman, J.R. Crutchfield and N.H. Packard, AppI. Phys. Lett. 37 (1980) 750. [2] D. D’Humieres, M.R. Beasley, B.H. Huberman and A. Libchaber, Phys. Rev. A 26 (1982) 3483. [3] A.H. McDonald and M. Plischke, Phys. Rev. B 27 (1983) 201. [4] W.C. Keer, M.B. Willians, A.R. Bishop, K. Fesser, P.S. Lomdahl and S.E. Trulliger, Z. Phys. B 59 (1985) 103. [5]W.J. Yeh and Y.H. Kao, Phys. Rev. Lett. 49(1982)1888. [6] E. Ben-Jacob, J. Goldhimsch, Y. Imryand S. Fishman, Phys. Rev. Lett. 49(1982) 1599. [7]M.R. Beasley and B.A. Huberman, Comm. Solid. State Phys. 10(1982)155. [8]M. Octavio, Phys. Rev. B 29(1984)1231. [9]M. lansiti, Qing Hu, R.M. Westervelt and M. Tinkham,
O.~,657
6.
—3.4
I
I
—3 2
—3
—2 8
I
I
—2.6
—2.4
—2.2
—2
Fig. 7. Log—log plot of the mean lifetime ofthe chaotic transient Tas a function of w,,—co.
observe a linear dependence with a slope ofthe order of —c5~—0.06, indicating a critical exponent 5ii~0.06.
6. Conclusions
Although in this paper we do not report on the calculations of Lyapunov exponents, power spectra analysis and phase portraits, we remark that they have all been considered and are in agreement with the types of behaviour that we have mentioned. From the results we have reported in this paper it is clear that the inclusion of the interference term produces substantial modifications ofthe dynamics. A more complete analysis could be done, both numerical and theoretical, for the interference model. Especially we mention that we are near to finish a study on the modifications suffered by the basins of attraction when varying the parameter y which controls the interference term from 1 to + 1, thus helping to investigate the change of the global behaviour of the system. —
342
Phys. Rev. Lett. 55 (1985) 747. [10] R.F. Miracky, M.H. Devoret and J. Clarke, Phys. Rev. A 31 (1985) 2509. [11] D.E. McCumber, Appl. Phys. 39 (1968) 3113. [12] A. Barone and G.J.Paterno, Physics and applications of the Josephson effect (Wiley, New York, 1982). [13]B.H. Tongue, Physica D 28 (1987) 401. [l4]Y.Yao,Phys.Lett.A 118 (1986) 59. [15] N.F. Pedersen, T.F. Finnegan and D.N. Langenberg, Phys. Rev.B6 (1972) 4151. [16]U.K.Poulsen,Phys.Lett.A41 (1972) 195. [17] D.N. Langenberg, Rev. Phys. Appl. 9 (1974) 35. [18]N.F. Pedersen, J. AppI. Phys. 44 (1973) 5120. [19] S.W. McDonald, C. Grebogi, E. Ott and J.A. Yorke, Physica D 17 (1985) 125. [20] C. Grebogi, E. Ott and J.A. Yorke, Physica D 7 (1983)181. [21] J.M.T. Thompson and H.B. Stewart, Nonlinear dynamics and chaos (Wiley, New York, 1986). [22] B.A. Huberman and J.P. Crutchfield, Phys. Rev. Lett. 43 (1979) 1743. [23] C. Grebogi, (1983) 935. E. Ott and J.A. Yorke, Phys. Rev. Lett. 50 [24] C. Grebogi, S.W. McDonald, E. Ott and J.A. Yorke, Phys. I,ett. A 99 (1983) 415. [25] C. Grebogi, E. Ott and J.A. Yorke, Physica D 24 (1987) 243. [26] E.G. Gwinn and R.M. Westervelt, Phys. Rev. Lett. 54 (1985)Thompson, 1613. [27] J.M.T. S.R. Bishop and L.M. Leung, Phys. Lett. A 121(1987) 116. [28] A.N. Lansbury and J.M.T. Thompson, Phys. Lett. A 150 (1990) 355. [29] J.A. Blackburn, H.J.T. Smith and D.E. Edmundson, Phys. Rev.A45 (1992) 593. [30] C. Grebogi, E. Ott and J.A. Yorke, Phys. Rev. Lett. 57 (1986) 1284.
PHYSICS LETTERS A
Crises and transient chaos in a driven Josephson junction F. Palmero Departamento de Fisica Aplicada, Universidad de Sevilla, P.O. Box 1065, Seville, Spain
and F. Romero Romero Departamento de Fjsjca Teorjca, Facultad deFIsica, Universidad deSevilla, P.O. Box 1065, Seville, Spain Received 26 January 1993; revised manuscript received 10 June 1993; accepted forpublication 14 June 1993 Communicatedby A.P. Fordy
Nonlinear dynamics of a single Josephson junction under radio-frequency excitation is analyzed by a model which takes into account interference effects and by means of numerical simulations. The study emphasizes the existence of a crisis, which is followedby chaotic transients over a portion of its state space. The mean lifetime is found to scale as I — I ~, where w is the frequency ofthe external excitation, a~,is its value at the crisis, and ö is the critical exponent.
1. Introduction Interest in the study of Josephson junctions has increased after the observation of Huberman et al. [1] that chaotic solutions might be responsible for the large noise increase in the Josephson parametric amplifiers [2]. Its mechanical analogue, the damped driven pendulum, has been the subject of extensive studies [3—5]and now it is clear that the Josephson junction is one of the powerful paradigms in nonlinear physics and that Josephson junction devices are imbued with a very rich dynamics and a variety of strange phenomena [6—10]. A typical Josephson junction can be described by McCumber’s model [11] and with a change of variables it is equivalent to the classical equation of motion of a simple pendulum with viscous damping, driven by a periodic external force. Nevertheless, to better describe the dynamics of the junction under the excitation by a radio-frequency current, the equation of motion must include a term that represents the interference effects between the quasiparticle-pair and quasiparticle currents, which is usually called the “cos 0” term. Most works about driven Josephson tunnelingjunctions have neglected
this interference term. The existing theories predict relative signs of the quasiparticle and quasiparticlepair currents which are in contradiction with some experimental results [12]. Then it appears reasonable to include in the evolution equations a “cos 0” term with a sign supported by experimental results. The analysis ofthe resulting equations might give rise to the prediction ofqualitatively new phenomena and might quantitatively affect the dynamical and fluctuation behaviour of Josephson devices. In this paper we have obtained computer solutions of the equation of motion modified by interference and analyzed the resulting data through various diagnostics tools: Poincaré sections, phase portraits, power spectra analysis, calculations of Lyapunov exponents and determination of basins of attraction. As the only varying control parameter we have considered the angular frequency of the external excitation. Thus, we have determined global bifurcation diagrams both for the noninterference model and for the interference model, which show substantial differences between them. By using the interpolated cell-mapping method (1CM) [131, we have obtained some basins of attraction, helping these to identify some topological
0375-9601/93/S 06.00 © 1993 Elsevier Science Publishers B.V. All rights reserved.
337
Volume 179, number 4,5
PHYSICS LETTERS A
types ofbifurcation. In this paper we are particularly interested in an interval of the control parameter where the system exhibits a chaotic blue sky catastrophe, associated with a hysteresis phenomenon, followed by transient chaos.
2. Model equation and numerical method McCumber’s model for a current-fed Josephson junction taking into account interference effects [141 has the following equivalent circuit equation, ~ + 4IEeR h (1 + y cos 0) 4xe dt2 0 = ‘o
sin (COrft)
+ I~sin 0
(1)
.
where 0 denotes the total phase difference across the junction, Cthe junction capacitor, h the Planck constant, e the electron charge, 4 the critical supercurrent, R0 the junction resistor corresponding to J~=ü, and y is the coherence coefficient which to 1o, of is therelated external the frequency, Wrt, and amplitude, excitation [15—171. This evolution equation can be transformed into the dimensionless form d2Ø ~
+a(l+ycosØ)~ +sinO=iesin(wr),
(2)
with T=(47~eI~/hC)”2t,a=R~’(h/4EeI~C)”2,~e= ‘o/j~~~ and w = w~/wo;wo = (4xeI~/hC)’/2~ In this form our equation has four control parameters, a, y, i~ and w, and can equivalently be written as
16 August 1993
considered values in the range 0.45 ~ w~0.73 and for ~e the fixed value ~e= 0.8. Some experiences on the Josephson plasma resonance in Pb—Pb tunnel junctions by Pedersen et al. [15] have shown a negative value for y, i.e. y= —0.9 ±0.2, then for our calculations we have chosen the value y= —0.8, which is in the range supported by the experimental results. Equations (3), (4) have been resolved in double precision arithmetic using a fourth-order variablestep-sizeRunge—Kutta method. In all cases, in order to eliminate transients, our solutions were run through over 200 periods of the driven current before actual data taking was started. For the numerical computations of the bifurcation diagrams the control parameter w has been increased (or decreased) in small steps. For each value of co we have obtained the steady-state solutions corresponding to a set of initial conditions randomly distributed in a determined subset of the phase space (X1, X2), and in a sufficient number such that all possible attractors are represented. Many aspects of the discussion are concerned with basins of attraction. Generally the determination of basins by numerical integration becomes extremely costly and time consuming. This approach has been used extensively in the literature when determining the basins of attraction of nonlinear systems [19]. Nevertheless, using the 1CM method one can determine the continuous time trajectory by spatially discretizing the system employing an interpolation procedure, with the advantage that one can generate an infinite number of distinct orbits from an Mx M mesh. In this paper we have used a 100 x 100 mesh over the area of phase space containing the attractors
dX 1 dr dX2 dr —
(3) =
—
a (1 + y cos X1 )X2
—
sin X1 + ~e sin (on), (4)
where the variable X2 is proportional to the difference of potential between the two superconductors. In order to obtain numerical solutions of eqs. (3), (4) we must choose a set of values for the parameters. Typical values of the coefficient a lie in the range 0.03 ~ a ~ 0.7, indicating that these systems are underdamped in the linear regime [11,18], and we have taken the fixed value a = 0.4. For w we have 338
of interest. For our model this grid is adequate to generate all the attractors that can be found by numerical integration, and the corresponding basins of attraction have enough resolution.
3. Bifurcation diagrams and crisis In order to appreciate to some extent the dynamical behaviour of the interference model as the parameter Co is varied, we have determined the global bifurcation diagram corresponding to the variable X2. Also, we have determined the global bifurcation diagrams for the corresponding noninterference model,
Volume 179, number 4,5
PHYSICS LETTERS A
16 August 1993
i.e. with the same values for all the parameters except for y, which now is 0. Figures 1 and 2 show both diagrams. From the comparison it is evident that the interference term modifies the Poincaré map rather substantially. Many aspects concerning the dynamics of the interference model could be analyzed, but
By decreasing the values of w from 0.73 we find that at w~0.672 a saddle node bifurcation occurs originating a stable period-three orbit and simultaneously an unstable orbit. This period-three solution is symmetric, i.e. it has inversion symmetry about the point (—It, 0) in the plane (X1, X2). Equation
in this paper we want to emphasize only some phenomena occurring in an interval around w=0.666.
(2) is invariant under the transformation X—~’ X, r—’r+ (2n+ 1)x/w, where n is an integer, then if X(r) is a solution of eq. (2), so is —X(t+ (2n+ 1)x/co). For the periodic solution encountered these two solutions are essentially the same, i e they differ by an integer number ofcomplete revolutions, so that
________________________________________ ~
.~
~
___,
—
(5)
—X(r+(2n+l)x/@)+2n=X(t)
As the period-three orbit is created there exist also chaotic solutions, these were originated after a cascade of period doubling bifurcations, i.e. a Feigenbaum route to chaos occurred for larger values of cv.
-
4
>~
~
•
..~.
~
‘~
~•‘.
At co= 0.666 a catastrophic event is about to appear. With extensive use of the 1CM method we have determined the basins of attraction corresponding to these coexisting attractors If we supenmpose the corresponding bidimensional Poincare sections and that ofthe unstable orbit, we obtain a composite por-
‘~ \
45
5
55
6
/
~ )-~
-.
65
7
Fig. 1. Global bifurcation diagram for the variable X2 correspondmg to the interference model with a=0 4 y= —0 8 i =0 8 andwintherange0.45~o~0.73.
.
~
‘
•~:~:‘
—2
—1
~
0
1
2
4’
3
Fig. 3. Basins of attraction in the plane (X1, X2) for the chaotic
I
.5
‘
xi
:
.45
~ ‘
,,
—3
•~frI
.
~
.55
I
,.,,
.6
.65
.7
w Fig. 2. Global bifurcation diagram for the variable X2 corresponding to the noninterference model, i.e. y=0, with a=0.4, —0.8, i,=0.8 and coin the range 0.45 ~ w~0.71.
solutions (blank) and the coexistent period-three solution (shaded) computed from eqs. (3), (4) for w=0.666, a=0.4, y=—0.8andi~=0.8.A1soshownarethePoincar~sectionsofthe chaotic attractor, the period-three orbit (centers of the full cirdes) and the unstable orbit (center of the open circle). The patterned-dot structure forced by the pixel separation ofthe printer must be ignored, since it is not caused by the numerics.
339
Volume 179, number 4,5
PHYSICS LETTERS A
16 August 1993
trait which can be observed in fig. 3. The basin boundary passes through the open circles, i.e. it belongs to the stable manifold of the unstable orbit, and we can see that it is about to collide with the chaotic attractor. In fact, at co= w,.y~0.66575, a tangency of the stable and the unstable manifold ofa period-three orbit, located on the basin boundary and accessible from the basin boundary ofthe chaotic attractor, occurs, originating a discontinuous disappearance of the chaotic attractor or crisis [20]. Thompson and Stewart [21] reported the existence of a similar phenomenon found in the forced oscillator equation d2X —0.25 dX —X+X3=A sin(w’r). (6)
terference model (eq. (2)), and in the range of cv where the period-three orbit coexists with the chaotic attractor, the basin boundary has also a fractallike structure, as can be seen in fig. 3. After the catastrophe, the basin of attraction of the chaotic attractor is incorporated to the basin of attraction of the period-three orbit, and when the symmetrybreakingbifurcation occurs the basin of attraction is divided in two regions, each one corresponding to a different orbit, yielding the structure shown in fig. 4. In the interpretation ofthis figure we must take into account that the pixel separation of the printer detennines a structure not caused by the numerics. Thus, the large ligthly shaded gray area is not an extension only of the densely shaded basin, but also includes “blank fragments”, and the resulting basin
This crisis is closely related to a hysteresis phenomenon observed by Huberman and Crutchfield [22]. For our case the hysteresis phenomenon is a necessary condition. When decreasing cv, destruction of the chaotic attractor is made possible by a previous creation of a saddle orbit, but as this is accompanied by the creation of a stable orbit, there exist simultaneously two attractors, the chaotic attractor and the period-three orbit. A symmetry-breaking bifurcation of the stable period-three orbit occurs for cv very close to w~,originating two different asymmetric orbits, each one being an inversion of the other with respect to the point (—x, 0) in the plane
boundary is also fractal. These facts determine the high sensitivity to initial conditions, not only when there exists a chaotic attractor, but when the only two coexistent attractors are periodic. Recent works relative to driven oscillators with cubic and quartic potential wells [27,28] showed some comparable forms of basin-breakup transitions.
-~—
I
I
_______
(X 1, X2).
.
4. Basins of attraction Another important question to consider is the one about the structure of the boundaries between basins of attraction when there is coexistence of attractors. Thus, Yorke and collaborators in a series of papers [19,23—25] treat the problem of how basins of attraction change as a system parameter varies continuously, and which are the qualitative changes of the basin boundaries as a parameter passes through certam critical values. They called such changes “metamorphoses”, and presented numerical evidences for fractal boundaries between basins of attraction for nonchaotic attractors. Fractal boundaries appear in many nonlinear systems as for example in the noninterference McCumber model [26]. For the in340
.-
________
_________________ _____
_____
.
.
:
-
Fig. 4. Basins ofattraction in the plane (X1, X2) for the two coexistent penod-three orbit for w = 0.664. The large lightly shaded area corresponds to the interspersion of the blank area and the densely shaded area to the basins ofattraction ofthe two coexistent orbits.
Volume 179, number 4,5
PHYSICS LETTERS A
5. Transient chaos
16 August 1993
~
...,. ~
Once the chaotic attractor has disappeared and for cv~w~,every trajectory starting from a given initial condition ultimately comes to one of the two coex istent period three attractors Nevertheless, the times employed to reach the final states differ and strongly depend on the initial conditions specially if they pertain to the basin of attraction of the destroyed chaotic attractor i e the system exhibits transient chaos This type of phenomenom has been observed in many different scenanos, as for example in the parametnc damped pendulum [29], where chaotic transients arise from the spatial congruence of the basinofattractionofastablestationarysolutionand the remnant of the destroyed chaotic attractor. Grebogi et al. [30] also showed some phenomena of transient chaos associated to bidimensional maps that appear both after a boundary crisis and after an interior crisis. For these cases the average lifetime Tof a chaotic transient scales as I p I where p denotes the control parameter and p~its critical value, and ô is the critical exponent ofthe chaotic transient. In order to find the distribution of the chaotic transient times for a given value of cv an appropiate set of initial conditions must be selected. As a first step we have taken 40000 initial conditions randomly distributed in the rectangle of the phase space _It~Xi~It,—2.5i~X2~l.5.By using the 1CM —
‘i~.,
~
I ..
,
‘~
...,,.
~
‘~.
~ .
—. “~
~
3
~J
~
~, ~.,
~.
-2
I I
1
.
0
1
2
3
Xi
Fig. 5. Union ofthe basin cores of the two coexistent attractors for w=0.665. ‘II’
—
— ~,
method we have followed their trajectories and selected all the initial points which reach a determined final state after a time of less than five periods ofthe external driving force. This selection gives rise to what is called “the basin core of the attractor”. For cv~0.665 there are two coexistent periodic attractors, the union of the two basin cores for cv = 0.665 can be observed in fig. 5 and as it can be seen is very similar to the basin of attraction ofthe periodic orbit previous to the chaotic catastrophe. To study transient chaos the appropiate set of initial conditions is determined by eliminating from the set of 40000 points all of them which pertain to the basin core. For each selected point we have determined the number of cycles required for the system to evolve to the final states. For each value of cv we have determined a histogram illustrating the number of initial conditions N versus the number of cycles of transient duration NT.
—
~
s
—
—
20
-
__________________________________ 40
60
80
100
120
N1
Fig. 6. Histogram of grouped data obtained for w= 0.665.
In fig. 6 we show the histogram corresponding to cv = 0.66 5. Inspection of this figure suggests that the probability distribution P( t) for the chaotic transient lifetime T is of the form P( ‘r) = (CI T) exp ( T), then T represents the mean lifetime of the chaotic transient. A plot of ln (N) versus ln (Ni) shows a straight line with slope T~. In order to determine the type of dependence of Ton cv in the proximity of cv,~,we have calculated the mean lifetime T for some different values of cv. Representing log10 ( T) versus log10(cv— cv,,), as it is shown in fig. 7, we can —
—
341
Volume 179, number 4,5
PHYSICS LETTERS A
16 August 1993
References
~
[1]B.A. Huberman, J.R. Crutchfield and N.H. Packard, AppI. Phys. Lett. 37 (1980) 750. [2] D. D’Humieres, M.R. Beasley, B.H. Huberman and A. Libchaber, Phys. Rev. A 26 (1982) 3483. [3] A.H. McDonald and M. Plischke, Phys. Rev. B 27 (1983) 201. [4] W.C. Keer, M.B. Willians, A.R. Bishop, K. Fesser, P.S. Lomdahl and S.E. Trulliger, Z. Phys. B 59 (1985) 103. [5]W.J. Yeh and Y.H. Kao, Phys. Rev. Lett. 49(1982)1888. [6] E. Ben-Jacob, J. Goldhimsch, Y. Imryand S. Fishman, Phys. Rev. Lett. 49(1982) 1599. [7]M.R. Beasley and B.A. Huberman, Comm. Solid. State Phys. 10(1982)155. [8]M. Octavio, Phys. Rev. B 29(1984)1231. [9]M. lansiti, Qing Hu, R.M. Westervelt and M. Tinkham,
O.~,657
6.
—3.4
I
I
—3 2
—3
—2 8
I
I
—2.6
—2.4
—2.2
—2
Fig. 7. Log—log plot of the mean lifetime ofthe chaotic transient Tas a function of w,,—co.
observe a linear dependence with a slope ofthe order of —c5~—0.06, indicating a critical exponent 5ii~0.06.
6. Conclusions
Although in this paper we do not report on the calculations of Lyapunov exponents, power spectra analysis and phase portraits, we remark that they have all been considered and are in agreement with the types of behaviour that we have mentioned. From the results we have reported in this paper it is clear that the inclusion of the interference term produces substantial modifications ofthe dynamics. A more complete analysis could be done, both numerical and theoretical, for the interference model. Especially we mention that we are near to finish a study on the modifications suffered by the basins of attraction when varying the parameter y which controls the interference term from 1 to + 1, thus helping to investigate the change of the global behaviour of the system. —
342
Phys. Rev. Lett. 55 (1985) 747. [10] R.F. Miracky, M.H. Devoret and J. Clarke, Phys. Rev. A 31 (1985) 2509. [11] D.E. McCumber, Appl. Phys. 39 (1968) 3113. [12] A. Barone and G.J.Paterno, Physics and applications of the Josephson effect (Wiley, New York, 1982). [13]B.H. Tongue, Physica D 28 (1987) 401. [l4]Y.Yao,Phys.Lett.A 118 (1986) 59. [15] N.F. Pedersen, T.F. Finnegan and D.N. Langenberg, Phys. Rev.B6 (1972) 4151. [16]U.K.Poulsen,Phys.Lett.A41 (1972) 195. [17] D.N. Langenberg, Rev. Phys. Appl. 9 (1974) 35. [18]N.F. Pedersen, J. AppI. Phys. 44 (1973) 5120. [19] S.W. McDonald, C. Grebogi, E. Ott and J.A. Yorke, Physica D 17 (1985) 125. [20] C. Grebogi, E. Ott and J.A. Yorke, Physica D 7 (1983)181. [21] J.M.T. Thompson and H.B. Stewart, Nonlinear dynamics and chaos (Wiley, New York, 1986). [22] B.A. Huberman and J.P. Crutchfield, Phys. Rev. Lett. 43 (1979) 1743. [23] C. Grebogi, (1983) 935. E. Ott and J.A. Yorke, Phys. Rev. Lett. 50 [24] C. Grebogi, S.W. McDonald, E. Ott and J.A. Yorke, Phys. I,ett. A 99 (1983) 415. [25] C. Grebogi, E. Ott and J.A. Yorke, Physica D 24 (1987) 243. [26] E.G. Gwinn and R.M. Westervelt, Phys. Rev. Lett. 54 (1985)Thompson, 1613. [27] J.M.T. S.R. Bishop and L.M. Leung, Phys. Lett. A 121(1987) 116. [28] A.N. Lansbury and J.M.T. Thompson, Phys. Lett. A 150 (1990) 355. [29] J.A. Blackburn, H.J.T. Smith and D.E. Edmundson, Phys. Rev.A45 (1992) 593. [30] C. Grebogi, E. Ott and J.A. Yorke, Phys. Rev. Lett. 57 (1986) 1284.