- Email: [email protected]
Prediction of bifurcations and chaos for an asymmetric elastic oscillator
Chaos, Printed
Sohtonr & Fractais in Great Britain
Vol. 2, No. 3, pp.303-321,
1992 0
0960.0779/92$5.00 + .cPJ 1992 Pergamon Press Lrd
Prediction of Bifurcations and Chaos for an Asymmetric Oscillator F. BENEDETTINI, Dipartimento
di Ingegneria
delle Strutture,
Elastic
G. REGA and A. SALVATORI
delle Acque
e de1 Terreno,
Universitl
dell’Aquila,
L’Aquila,
Italy
Abstract-The nonlinear chaotic response of a harmonically forced elastic oscillator with quadratic and cubic nonlinearities is studied. The possibility of obtaining reliable prediction of bifurcations and chaos through combined use of stability analysis of low-order approximate solutions and accurate localized point-by-point and cell mapping computer simulations is examined. Some satisfactory results are obtained for the bifurcation predictive capability and the location of regions where chaos actually occurs.
1. BACKGROUND
AND MOTIVATION
Oscillators with asymmetric and symmetric nonlinearities have received increasing attention in recent years [l-7]. Besides the theoretical interest associated with the presence of both kinds of nonlinear terms, they have been specifically referred to in structural applications as representing the dynamic behaviour of elastic systems with initial curvature undergoing finite displacements. Systems governed by an equation of motion with quadratic and cubic nonlinearities in the form 4 + /Uj + c1q + c2q2 + c,q3 = PcosQt
(1)
have most often been studied. In Refs [4, 6, 71 a dimensionless oscillator of type (1) was considered as a one degree-of-freedom model to describe the planar finite dynamics of a prestressed suspended cable subjected to vertical harmonic forcing. Its coefficients (cl = 1, c2 = 35.953, c3 = 534.53) were given values representative of the properties of a cable used in overhead transmission lines and characterize it as a strongly nonlinear oscillator with a single equilibrium position. In particular, it was shown how, notwithstanding the model’s simplicity, the presence of even and odd nonlinear terms gives rise to a rich and complex set of responses [7]. A response chart in the forcing parameter space (Q, P) was built through point-by-point time integration of the motion equation with fixed initial conditions (Fig. 1). Various periodic solutions occur in the three main zones into which that space can be subdivided, namely the neighbourhoods of the l/2- and l/3-subharmonic resonance conditions (associated with the two kinds of nonlinearities, respectively), and that of the superharmonic resonance conditions. In each zone, with proper values of the forcing amplitude, the periodic solutions bifurcate to chaos following various routes and different kinds of strange attractors occur. The actual chaotic nature of those responses was tested using various qualitative and quantitative dynamic measures [7]. From the theoretical point of view, there is strong interest in obtaining reliable predictions of the conditions beyond which the periodic orbits are destroyed and persistent, steady-state strange attractors are established [8-lo]: unfortunately, such an objective is hard to obtain, so practically all developed criteria lead just to necessary conditions for 303
chaos. This is connected with the strong sensitivity of the response to control parameter variations and with the associated fractal nature of the basin boundaries not only in the initial condition space but also in the system parameter space. Thus, it is always possible to find subsets of parameter values leading to responses which look like irregular ones but cannot be directly classified as chaotic. However, such a circumstance is not a great drawback from the point of view of engineering system applications, for which the simple assessment of unpredictability of the response-or the onset of a dynamic regime different from the ordinary one-is often more important than its precise characterization as chaotic. With these problems in mind, it seemed worthwhile to examine the possibility of predicting the disappearance of the periodic motions obtained for system (1) and of locating limited regions in the system parameter space where irregular motions and chaos can occur, while concurrently analyzing the effect of variations of the i.c. Among the various predictive criteria-either based on the homoclinic orbit bifurcation [S-lo] or of heuristic type [9, 11,121 or associated with the stability analysis of regular nonlinear oscillations [l, 3,5,6, 13-16]-the latter were selected according to suggestions first made by Szemplinska-Stupnicka [l]. Indeed, they are of special interest in the present case because they allow the prediction of bifurcations which are possible precursors to chaos. and the establishment of relationships between the irregular and the dominant regular solutions occurring in various resonance zones. With these criteria, approximate steady-state solutions for the system are first obtained through one of various analytical methods and then the stability of the solutions is evaluated. One major procedural problem lies in the choice of the type and number of harmonic terms to be used for obtaining both accurate solutions and stability analysis. Approaches involving low [l, 3.5,6,13] or high-order [14-161 approximations can be distinguished, the latter are able, in principle, to furnish better results. However, two points are considered: (i) the main advantage of developing and applying criteria of this type is the possibility of recovering the huge patrimony of knowledge developed in the classical field of nonlinear oscillations for completing the understanding of new dynamic phenomena; (ii) interest is focused on identifying zones of response unpredictability rather than strictly chaotic motions. Within this ambit, it is worth examining how much accuracy can be obtained in the prediction of irregular responses by means of analytical solutions that are as simple as possible. Thus, low-order approximations will be considered later. Use of higher approximations is likely to be necessary in the study of more complex systems or phenomena, mostly as far as reliable stability evaluations are concerned [16], but these can sometimes be useless for engineering purposes due to roughness or uncertainty of the model considered. In the next section, the bifurcation predictive capability of approximate solutions properly chosen in the various regions of resonance is analyzed with the help of the results of accurate but localized computer simulations. 2. APPROXIMATE
PERIODIC SOLUTIONS, STABILITY ANALYSIS, BIFURCATIONS, OF ATTRACTION
BASINS
Equation (1) being asymmetric, its period jT solution q(t) does not satisfy the symmetry property, i.e. it is q( t + jT/2) # -q(t), j = 1, 2, 3, . . . . If one seeks a solution represented by finite Fourier series, an asymmetric solution generally contains a constant term and both odd and even harmonic components:
Prediction
of chaos for an asymmetric
305
oscillator
(2) can be determined through, e.g. the harmonic The coefficients of/j, ~i~j, in equation balance method. Examination of local stability of the steady-state period jT solution is made by adding to it a small disturbance, g(t) = q(t) + 6q(t). Then, inserting it into equation (1) and neglecting nonlinear terms in 6q(t), a linearized Hill’s type variational equation is obtained, 64 + ,LL& + g(q)dq
= 0,
(3)
in which the periodic coefficient g(q), evaluated through equation (2), has period jT. The Floquet theory [17] calls for solutions to this equation of the form (Floquet multiplier) of the monodromy 6q(t + jT) = A6q(t), where h. is an eigenvalue matrix C associated to a fundamental matrix solution Q, of equation (3), @(t + jT) = C Q(t). By assuming Q(O) = I, it results C = @(jT). The matrix Q(jT) can be computed by numerically integrating equation (3) for fixed values of P and S2 in the interval [0, jT] subject to two sets of i.c., (i) 6q, (0) = 1, 64, (0) = 0 and (ii) 6q, (0) = 0, 84, (0) = 1. The values of i\. determine the stability of the approximate solution q(t). IAl< 1, which means eigenvalues of C remaining inside the unit circle in the complex plane, is a necessary When varying a control parameter of the system, the and sufficient stability condition. manner in which d leaves the unit circle characterizes the local qualitative bifurcations + 1, one has occurring to the periodic orbit. When it leaves the unit circle through 6q(t + jT) = 6q( t). Thus, q(t) is a jT-periodic attractor and the system undergoes a jump bifurcation. On the other hand, when A leaves the unit circle through -1, it results 6q(t + 2jT) = as(t). Thus q(t) is a 2jT-periodic attractor and the system undergoes a period-doubling bifurcation.
2.1.
Neighbourhood
of l/2-subharmonic
resonance
The computer simulations of equation (1) with zero i.c. made to obtain Fig. 1, showed that the route to chaos from the right in the neighbourhood of the 1/2_subharmonic resonance is characterized by a dominant sequence of p.d. bifurcations [4,7]. This is clearly observed in Fig. 2, where the period chart of the response obtained at a fixed value of forcing amplitude (P = 0.04) is reported. To examine this phenomenon within the framework of stability analysis of approximate solutions, the underlying dominant period 1 response is described by the simplest asymmetric expansion, containing just the drift and fundamental harmonic: q(t) = a0 + alfcos(Q2f Substituting it into equation (1) and applying system of three algebraic nonlinear equations: 3c,(a2,alf
+ afJ4)
+ 8).
the harmonic
+ 2czaoalf + alf(l
(4) balance
+ 3a,a$)
leads
to a
- Q2) = Pcos6
alfpS2 = - Psin6 c,(2ai
method
(5)
+ c2(2ai + aff) + 2a. = 0.
By solving it using a numerical procedure, the unknown amplitudes and phases of the solution are obtained; in particular, the frequency-response curve can be built. Analysis of its local stability is made following the procedure outlined previously. The eigenvalues of C have been repeatedly computed as Q is varied at fixed values of P. The frequency-reszones of instability (thin lines) are ponse curve for the amplitude aIf and the relevant plotted in Fig. 3(a) (P = 0.04). Besides the classical unstable hysteresis branch, a further
306
F. BENEDETTINI et al, 1.0
0.0
1’
0.4
0.4
0.3
0.3
0.2
0.2
0.1
0.1
0.0
---1
1.0
0.0
I-3
Period
Period
1 m
6xN
Fig. 2. Period
chart
1.70
I
3.0
Period
m__
Fig. 1. Periodic
1.66
I
2.0
1.74
2 i
Period
and chaotic
I .7a
in the neighbourhood
response
1.82
3
in control
1.86
of l/2-subharmonic
R
Period
m
Transition
0.0
I
4.0 6
Chaos space (P, n).
1.90
resonance
1.94
1.98
(P = 0.04, zero i.c.).
Prediction
of chaos for an asymmetric
oscillator
307
large zone of instability is attained through p.d. bifurcations on accordingly, the establishment of a period 2 solution is expected. Thus one can consider the simplest expansion of period 2, namely q(t) and obtain amplitudes
= a0 + a,cos(Qr
+ 79) + a,,,cos(Q/2~t nonlinear method,
,/+os(6 ($3% + $4Q2
+ 2C,U” + &z:
- a,(@ (ic3u0 (3c,u,u,
+ ic2)a&sin(t9
+ c2ul)ul/2COS(19
-
lower
branch:
+ q),
the following system of five algebraic and phases through the harmonic balance - 299 + a,(3c3a;
the
equations
(6) in the
unknown
+ ;C&)
1) = P cos 6
- 2cp) - a,,~&? = Psinff
- 2q) + 2c*a~~~,,/ +
c3
i
3u; + ; u: + ; a:,2 a,/2 31
!a2 + 1-qu:/2=0
1
i (3c3uoul ~c3u&c0s(6
+ c~u~)u~,~ sin (79 - 2y) -
+LZ~/~@ = 0
- 2l7I) + c*[ui + ;
and solve it numerically. At the considered value of forcing amplitude, the range of existence of this solution coincides precisely with the range of second-order instability of the period 1 solution: the frequency-response curve for the amplitude ali of the subharmonic component is plotted in Fig. 3(a), as well. By repeating the stability analysis, the lower branch is seen to be practically always stable, while the upper branch shows a large zone of instability attained again through p.d. bifurcations: the occurrence of a period 4 solution is thus expected.
(4
(b)
cl.0
0.5
1.0
Forcing Fig. 3. Frequency-response
15
2.0
2.5
3.0
3.5
4.0
frequency curves: approximate solution (thick, stable and thin, unstable) lation results with own i.c. (black, also zero i.c.). (P = 0.04).
and computer
simu-
308
F. BEN~D~TINI et al
Use of higher subharmonic approximations would make likely the possibility of further p.d. bifurcations evident. However, these occur in very small ranges and the desire is to limit the analysis of the bifurcation predictive capability to low-order approximations. So, the predictions obtained are compared directly with the results of computer simulations. The amplitudes of the harmonic and 1/2subharmonic components obtained by calculating the FFT of the simulation results are compared in Fig. 3(b) with the corresponding harmonic balance curves. Different from the stability analysis, numerical simulations require specification of the i.c. The values q(O), Q(O), obtained by inserting the harmonic balance amplitudes and phases in equations (4) and (6), are considered here (white and black dots in Fig. 3(b)) compared to the zero values already ascertained for obtaining the chart in Fig. 1 (black dots in Fig. 3(b)). It is observed that by using the conditions defined by equation (4), the upper a,, branch is almost completely matched, while its lower branch is matched only where it is stable according to the approximate analysis. This is not shown in Fig. 3(b). Where the analysis indicates occurrence of p.d. bifurcations, the conditions allow the curves a,‘?, a, of the period 2 solution to be followed. When using i.c. defined by equation (6) and decreasing the frequency, the lower u,/? branch is matched from the frequency value (Q G 2.26) where the period 2 solution arises up to a point close to the vertical tangent. This is consistent with the branch’s observed stability and, thus, agreement with the theoretical curve is very good. The response then jumps to the upper ali branch. The upper branch is followed from its peak amplitude up to just the point where the stability analysis shows the occurrence of a p.d. bifurcation. Left of this point, a l/4-subharmonic component appears in the numerical solution and then successive p.d. bifurcations to period 8, 16, . . are observed. As the bifurcation proceeds towards chaos, the amplitudes of 1/2_subharmonic component decrease slightly with respect to the theoretical values, due to the presence of further subharmonic components. The zone of chaos extends approximately between Q = 1.82 and Q = 1.70. Below this value, even with the i.c. defined by equation (6), the response jumps to upper period 1 and the relevant ~7,~ branch is followed up to a frequency value (Q z 1.52) where the response jumps down to period 2 and the corresponding theoretical curves-which are now stable again-can be followed numerically. The analysis confirms substantially the points made in [I], where the same problem was studied for a system with asymmetric nonlinearity similar to equation (1) by considering different values of the system parameters. Slight differences are represented by the occurrence of the chaos zone between the two points of vertical tangencies and by the possibility of moving along the al:? branch left of that zone: the latter suggests that the irregular motions found are actually associated with the specific second-order instability of the period 2 solution. When decreasing the frequency with zero i.c., the lower a1;1 branch is followed up to a value of frequency (Q = 1.96) below which the process of bifurcation to chaos develops and solutions with higher period (e.g. period 4 at B = 1.92) are established exhibiting the 1/2subharmonic component with upper amplitude. In substance, the jump to the upper branch has shifted to a slightly higher value of frequency with respect to the previous case. These local shifts heuristically suggest the possibility of a fractal-like nature of the boundaries between different periodic (and/or chaotic) regions in the forcing parameter space depending on the i.c. The capability of stability numerical analysis of the approximate solutions to accurately predict bifurcations and delimitate regions where chaos can occur is illustrated more limit curves of generally in Fig. 4 in the forcing parameter space. Here, the stability approximate period 1 (equation 4) and period 2 (equation 6) solutions and the locus of
Prediction
of chaos for an asymmetric
309
oscillator
points of vertical tangency to the period 2 frequency-response curve are superimposed on the regions of response obtained through numerical simulations around the l@-subharmanic resonance. Of course, for a correct comparison of numerical and theoretical results, one must remember that the former refer to zero i.c. These belong in various zones: to the domain of attraction of the lower period 1 solution for frequencies higher than the corresponding right stability limit curve (zone 1 in Fig. 4); left of this curve, to the domain of lower period 2, that is always stable (zone 2); where lower period 2 does not exist, to the domain of upper one (zone 3); finally, to the domain of upper period 1 which is always stable (zone 4). Embedded in these zones and even nearly bounded on the right by the limit curves of first or second order instability of period 2 solution, a large zone of irregular (highly periodic) and chaotic motions occurs. To clarify all features of the system dynamic behaviour in the considered range of forcing parameters, one has to complement the regions of periodic and chaotic response of Fig. 4 with systematic construction of basins of attraction in the i.c. space. In particular, the extension of the domain of chaotic solutions must be checked. Shown in Fig. 5(a-f) are the basins of the different solutions obtained with a cell mapping technique [18] for various values of decreasing frequency at the same forcing amplitude as in Fig. 3. To the right of the zone of existence of the period 2 solution, the basins of the upper (light blue) and lower (blue) period 1 spiral into each other. At the value Q = 2.26 (a), corresponding to the peak amplitude of the ali curve, very thin tongues (yellow) relevant to the upper period 2 solution appear within the basin of lower period 1. Just to the left of that value, consistent with the p.d. bifurcation prediction of this latter solution, the blue basin corresponds to lower period 2, though preserving the same knob shape. Moving towards lower frequency values, the two basins of lower and upper period 2 solutions evolve differently: the former one vanishes while the upper one increases, though
l .
.
l
Stability
limit
of
period
1
solutron
o
Stability
limit
of
period
2
solution
-
Limit
of
vertical
tangency
in
sub.
l/2
i 0.1
i
Period
El
, Ij ‘$
----I I!zrI
Transition
Penod
1.3 Fig. 4. Bifurcation
2.0 predictive
capability Comparison
1
to
3
Period
i
2
chaos
Chaos
2.6 of stability numerical analysis of period 1 and 2 approximate with computer simulation results (zero i.c.).
solutions.
310
F. BENEDETTINI et al
corresponding to solutions with higher and higher periodicity which eventually end in chaos (b-d). The strong fractal nature of the relevant basin boundaries is worth observing (b, c). At the same time, the size of the knob diminishes and the basin of attraction of upper period 1 solution increases notably. The ‘sudden’ change of the response from chaotic to upper period 1 occurring at the left boundary of the chaos zone in Figs 3(b) and 4. corresponds to strong mixing of the relevant basins in the i.c. space, that however is observable only through use of a much narrower window (e). To the left of that zone, any i.c. leads to the upper period 1 solution up to the frequency value (Q z 1.52) where a new spiraling tongue, corresponding to solutions with higher periodicity, appears in the i.c. space and the period 2 solution is again established through a reverse p.d. bifurcation (f). As the frequency is further decreased, this tongue becomes thinner and thinner and when the period 2 solution disappears, it corresponds to the lower period 1 solution still occurring in a short frequency range (see Fig. 3(b)). The important role played by the i.c. clearly emerges. Moreover, it is worth noting how the observed instability of the dominant approximate periodic solution-though a necessary condition for chaos-also suggests the most likely i.c. to consider for actually obtaining it. In the present case, they were rather successful.
2.2.
Neighbourhood
of l/3-subharmonic
resonance
subharmonic and chaotic responses According to the results of computer simulations, arise in this frequency range following more mingled routes with respect to those occurring in the neighbourhood of l/2-subharmonic resonance. This can be observed in the period arise at certain frequency chart of Fig. 6, again obtained with zero i.c. Period 3 solutions values and intertwine with period 6 solutions into which they bifurcate themselves. The latter give rise to several responses with period multiple of 6 increasing irregularly up to very high values. Accordingly, the chaotic strange attractor exhibits six independent bundles orginating from the period 6 solution around which the motion fluctuates [7]. Besides the fundamental solution equation (4), approximate steady-state solutions are sought by considering the simplest expansion of period 3, namely q(t)
= a,, + a,cos(Qt
+ 79) + a,/,cos(Q/3-
f + v).
_.~~_
1
Ch.
4
Fig. 6. Period
chart
in the neighbourhood
of l/3-subharmonic
resonance
(P = 0.40, zero i.c.)
(8)
Prediction
of chaos for an asymmetric
oscillator
Fig. 5. Basins of attraction in the neighbourhood of l/2-subharmonic resonance. P = 0.04; P = 2.26 (a), 2.00 (b), 1.92 (c), 1.76 (d), 1.69736 (e), 1.50 (f). Range for q(O), Q(O): -0.2 + 0.2 (-0.02 + 0.02 in case (e)).
312
Fig. 8. Basins
F. BENEDE~INI
etal
of attraction in the neighbourhood of l/3-subharmonic resonance. P = 0.40; B = 3.35 (a), 3.42 (b), -0.3 + 0.3; for Q(0): -0.4 + 0.4. 3.58 (c), 3.82 (d), 4.13 (e), 5.00 (f). Range for q(0):
Prediction
of chaos
for an asymmetric
oscillator
313
0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2 0.08
0.08
0.07
0.07
0.06 g g 0.05
;..
&.. . 8004 ! _4 zi g ................_./_..._....._ .. ..
I . . . . . . .
...)
eo.03
‘_ j
.
.......
0.06
\ ......
0.05
(
......
0.04
,
......
0.03
2
0.02
.
.
. .
: ......
.
0.01
0.02
0.01
0.00 0.00 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.3 0.3 1.0 1.1 1.2
"A,
EA,
Fig. 11. Frequency-response
Xi-A,OA, A,*A6.A,~A,AA,
Forcing frequency
curves: approximate solution (thick, stable and simulation results with own i.c. (P = 0.10).
thin,
unstable)
and
computer
314
F. BENEDE~INI
et al.
Pig. 13. Basins of attraction in the neighbourhood of order 2 superharmonic (b), 0.89 (c), 0.90 (d), 0.91 (e), 0.93 (f), 0.954 (g), 0.95526 (h). Tl green,
peak. P = 0.10; Q = 0.85 (a), 0.87 T2 yellow, T3 light green, T4 blue.
Prediction
of chaos for an asymmetric
By applying the harmonic balance method, in the amplitudes and phases is obtained 14 c&c0s(6
- 37#) + a,(3c,a;
2C$z~ + +ju:
+
i c,aT,3sin(7? 3 4
the system
315
oscillator
of five algebraic + +&j)
nonlinear
- a@
equations
- 1) = Pcos6
- 3rjJ) - alpQ = Psin6
3 - 3r/~) + 2c2a0al/3 + cx 3ai + 5 a: + 43 ali3 2)al,~+(l-~)&?l,3=0 (
c,alf&cos(6
: c3ala:j3sin c2[ai + ;(a:
(6 - 3q)
+ f a,ppQ
= 0
(9)
+ ~75~)] + c3[ai + tao(a: + ~z:,~)] + a0 = 0,
and then solved numerically. The frequency-response curves for the period 1 solution and for the amplitude ali3 of the period 3 solution are plotted in Fig. 7(a); a value of forcing amplitude (P = 0.40) is considered for which the simulations showed occurrence of chaos (Fig. 6). The stability analyses show that both the upper and the lower branch of period 1 are always stable in the range of interest, while the upper branch of the subharmonic component exhibits a large zone of instability (3.65 < Q < 4.98) attained through p.d. bifurcations: at the branching points, the occurrence of a period 6 solution is thus expected. Again, the predictions are compared directly with the results of computer simulations obtained with both zero i.c. (black dots) and the conditions corresponding to the harmonic balance approximations (white and black dots) in Fig. 7(b). With those defined by equation (4), the two stable branches of aIf are completely matched except in a very small part of the lower one near the vertical tangency to the aI/3 curve. With the conditions defined by equation (S), period 3 responses are obtained with very high accuracy in the whole range of theoretical existence of l/3-subharmonic solution, except in a narrow zone (3.60 < Q < 3.74, the same found when using zero i.c., see Fig. 6) where period 6 solutions appear and, correspondingly, the amplitude a,/3 decreases slightly with respect to the theoretical value. The responses obtained with zero i.c., already shown in Fig. 6, are also
(a)
(b)
0.0
2.0
Forcing Fig. 7. Frequency-response
40
6.0
frequency
a.0
10.0
IL
3.00
~~~~~
~
4.00
Forcing
5.00
6 00
frrquency
curves: approximate solution (thick, stable and thin, unstable) lation results with own i.c. (black, also zero i.c.).(P = 0.40).
and computer
simu-
F.BENLD~TTINI et al
316
indicated from the black dots and the arrow in Fig. 7(b), where the frequency range in which chaos occurs is also shown. The bifurcation from period 3 to period 6 response evident from previous results is consistent with the bifurcation prediction of the stability analysis of the period 3 approximation (8). However, the notably different kind of response (period 3 or chaotic) obtained with different i.c. raises some doubt about the overall reliability of the analysis and suggests a strong influence of the i.c. on the results. Better understanding of the complex set of responses which occur can be achieved by examining the evolution of the relevant basins of attraction with increasing frequency (Fig. 8(a-f)). To the left of the zone under investigation, only a small basin of lower period 1 (violet) occurs within upper period 1 (red) (a). Where the stable period 3 solution arises, thin tongues (yellow) enter the former basin (b), and then quickly expand over practically the whole basin. Thereafter, consistent with the prediction of p.d. bifurcation of the relevant stability analysis, this basin corresponds to solutions with increasing period (6, 12, . .) ending in chaos (c, d). Two points are worth noticing. First, the basin of attraction of lower period 1 is reduced up to the knob core, but it never vanishes: this is consistent with the overall stability of the solution shown by the approximate analysis. Second, the complex evolution to chaos of the response within the knob is accompanied by the onset of new tongues of period 3 solution (pink) surrounding the core of lower period 1 basin (d, e). Thus, the period 3 response obtained for its own i.c. within the chaos zone and charted in Fig. 7(b), is associated only with the existence of these tongues, to which those conditions belong. However, though taking into account the strong mixing of the various basins occurring within the knob, most of the relevant i.c. lead to responses with period other than 3 (lower periods 1, 2, 4 and 6 in Fig. 8(e)), thus confirming in general terms the instability prediction of that solution obtained with the approximate analysis. Such a circumstance is understood more clearly through comparison of the detailed picture of the response regions obtained with computer simulations in the zone of main interest around the l/3-subharmonic resonance, and the regions of stable and unstable period 3 solutions given by the approximate analysis (Fig. 9). A shift of the approximate unstable region towards higher frequencies with respect to the numerical one is observed. But the overall picture is very similar in the two cases, with the wedge of stable solutions entering the upper part of the figure. Moreover, the shift has limited meaning since the numerical results correspond only to zero i.c. On the other hand, the basins of attraction obtained at different forcing amplitudes show the period 3 solution tends to occupy the whole knob in the upper part of Fig. 9, consistent with its approximate stability. Instead. it is greatly reduced by the occurrence of various further solutions-with marked fractal character of the knob itself-in the lower part of Fig. 9, where the period 3 solution practically vanishes (e.g. it does not exist at P = 0.2, Q = 3.25). Finally, when considering the even greater role played by the i.c. in this parameter zone, the capability of the simplest approximate solution to delimitate regions where complex and chaotic motions can occur is satisfactorily established. Coming back to Fig. 8, it is observed that with the further increase of the frequency (f), the basin of period 3 solution (yellow), now stable again, increases notably and remains the only one in the knob-together with the lower 1 basin-while the knob itself tends to occupy the whole window.
2.3.
Neighbourhood
of superharmonic
resonances
the response picture is quite rich and complex with To the left of primary resonance, many periodic and chaotic zones not so well defined, but following a certain regular pattern
Prediction
of chaos for an asymmetric
317
oscillator
. .
0.46
.
1 .
.
limit 3
Limit of vertical tangency in sub
-
.
Stability of period solution
. .
0.36
3 .3 Fig. 9. Regions
3.5
3.7
of stable
3.9
4.3
4.1
and unstable
4.5
4.9
4.7
5.1
o
Period
1
o
Period
3
[1
Period
1,3
D
Period
1,3,6
u
Period
6
o
Period
6*N
-
Chaos
5.3
period 3 solutions with approximate simulation results (zero i.c.).
analysis.
Comparison
with computer
of inclined and nearly parallel stripes embedded within the dominant period 1 solution and extending over the entire range of forcing amplitude (Fig. 1). Several transition zones occur in which the response periodicity is notably sensitive to small variations of forcing parameters (see also the period chart of Fig. 10, still obtained with zero i.c.). It is thus easy to predict high sensitivity to variations of i.c., too. An approximate steady-state solution is sought by considering a period 1 solution improved by the presence of higher harmonics, namely the order 2 and 3 superharmonic terms: + p) + a,cos(3Qt
q(r) = a(J + UlCOS(Qf + ?Y)+ azcos(251t
__
Chaos
_._
+
I@).
_--I-I_
-
72 a
0
‘f: g
20
Qi z k !.G p:
aL 4 1
-_A
0.2 Fig. 10. Period
.:--y.:---: .. . .
chart
0.3
.. .
__--
_. ___
0.4
in the neighbourhood
.0.5
‘t_
0.6
of superharmonic
_
-
R
-
0.7 resonances
---~ 0.8
0.9
(P = 0.10, zero i.c.).
(10)
318
F. BENEDETT.INIet al
With the usual technique, 3c+ 4
[afcos(38
the system
- I/J) + a:cos(6
of seven algebraic
nonlinear
equations
- 2p + q)] + az(cz + 3c,a,)[a,cos(26
reads:
- cp)
3c3 a: + fp - q)] + a, ao(2c2 + 3cjao) + y y + ai + ai i [ I
+ a3cos(6
1 43c3a3
= Pcos-9
+ (1 - n*)
[afsin(38
-
+ a,sin(@ al(c2 + 3c,a,)
q) + azsin (6 - 2~ + q)] + a2(c2 + 3c,ao)[aI sin (2tY - q) q)] - ,@a,
+ q -
- q) + a3 cos(B + cp - tp> + 3c,a,
$cos(211 I
1 al(c2 +3c3ad [ + ai + 2
?
+ 2
7-
ai + 4
Tsin
-
3c3ala2a3
2
+ a,(1 - 4Q2 + 2c2a,)
(21) -
__
2q1-t
l/Q)
= 0
q) + a3 sin (8 + q -
1
(11)
ye)
sin(8-2fp++)-2yQa2=0
3
7cos(31’)
ala3 1 I2cos(6-
= Psint9
- v) + ala2 (c2 + 3c,a,,)cos(6 i
+ i 22 3c,a3 ai + :
+ :
3c3a7 + y - I/I) + -ycos(6
- 2y, + q)
ai + 4 i + a,(1 - 9R2 + 2c,a,)
1
= 0
3 ysin(3B
3c3ala2
-
ijj) + ala2 [ (c2 + 3c,a0)sin(6
+ 3pQa3 = 0 1 ~cos(213-cj7)+~cos(r’l+~-~) + a0(c2aC, + c,ai
+ q -
1
v) - Tsm(8 3c3a2
i +(cZ+3c3a,,)~+$+$
2
- 295 + I./I)I 2
2
1
+ 1) = 0.
The frequency-response curves obtained for a,, a2, a3 at the same value of forcing amplitude as in Fig. 10 are plotted in Fig. 11. The considerable modification of curve a, of the fundamental harmonic associated with the two superharmonic resonances Q = l/2 and Q = l/3 is clearly observable. Correspondingly, the superharmonic amplitudes a2 and a3 increase strongly [19,20]. The construction of the curves through numerical solution of system (11) was not forced upward to obtain the unstable branches of the loops which are likely to occur between the points of vertical tangency for both the harmonic and the superharmonic components, the former downward and the latter upward [19], but they are strongly suggested by the trends of the curves themselves. The stable parts of the approximate solution are denoted as usual with thick lines in Fig. 11, but these are limited only to the a, curve for the sake of clarity. In the same figure the results obtained from computer simulations, made on a very fine frequency mesh by
319
Prediction of chaos for an asymmetric oscillator
considering i.c. corresponding to harmonic balance approximation (lo), are also plotted. The overall agreement is satisfactory. Discrepancies occur: (i) as far as the solution is concerned, mostly in the lowest frequency range, where higher superharmonics (namely, 4 and 5)-obtained from computer simulations and also plotted in Fig. 11-play an important role; (ii) as far as stability is concerned, with respect to the localization of some unstable zones within the frequency range. Consistent with the findings from Fig. 10, the picture of motions ensuing from the basins of attraction is very complex. Thus it is not possible to follow their evolution, but only to focus on some local results. Three main frequency zones can be distinguished in Fig. 11, two located in the neighbourhoods of the order 3 and 2 superharmonic peaks, the third one just between them. As to the first zone, it can be seen that to the left of the peak (Q G 0.52), the i.c. corresponding to both upper and lower approximate period 1 solution belong to the upper 1 basin: this is consistent with the stability and instability predictions obtained respectively for the two solutions (Fig. 11). By slightly increasing the frequency, a chaotic basin-to 1 belong-spirals within the upper period 1 basin which the i.c. of lower period (Fig. 12(a)), that decreases (b) till the solution itself vanishes (vertical tangency at Q G 0.5805 in Fig. 11). At the same time, the solution in the growing basin-which eventually remains alone-evolves quickly and irregularly through many different degrees of periodicity. In the second zone, corresponding to the intermediate branch of fundamental harmonic in Fig. 11, a rather variable response picture occurs in the i.c. space, with one or more basins of low periodicity. Between about Q = 0.7 and 0.8, corresponding to a decrease of basin fills in the overall i.c. window. the simulation amplitude al, a chaotic or subharmonic It is worth observing how the approximate analysis reveals a zone of unstable period 1 solutions, but this is shifted towards higher frequency values with respect to the simulation results. Accounting for higher superharmonics in the approximate solution would be likely to be necessary for better stability prediction [16]. The third zone corresponds to the order 2 superharmonic peak and is quite interesting. The basins of attraction with the increasing of frequency (Fig. 13) evolve from a single
(b)
Fig. 12. Basins
of attraction
in the neighbourhood (b). Range
of order 3 superharmonic for q(O), 4(O): -0.2 + 0.2.
peak.
P = 0.10;
Q = 0.56 (a), 0.58
period I basin (?‘I) to four coexisting hasins (;I-d)--- all of period I --appearing in (urn (‘1’2. ‘1’3, ‘1’4). For SZ ~~0.90. the nunihcr of period I solutions diminishes down to one’ (c-h) hut with ;I vanishing scqucncc different from the appcarancc scqucncc. Indeed, the surviving basin corresponds to ;I motion (‘1‘3) other than the original enc. ‘I’hcsc two solutions must Iw rclwtcd to the left and right hranchcs of the I’untla~nc~~t~~l in f;ig. I I. A clear undcrstanclil~g of what happens is xhicvctl
h;~nnonic through
(I, around the peak a comparison ol’ I:igs
I3 and 14. The latter shows an cnlargcmcnt of the 0, curve around the peak, ohtainctl from computer simulations with proper i.c. 13esidcs the two main period I branches--‘I’1 from the left and T3 to the right- which overlap and arc likely to co~~ncct with WC‘ another below via an unstahlc loop 1201, two more period I brnnchcs (‘I’2 and ‘1‘4) arise in ;I small frequency range and exhibit nearly the same (7, amplitutlcs as the previous ones. ‘Thcsc two branches overlap in a short range as well. where four stahlr solutions exist, and at-c’likely to connect with one another via unstable branches that could he dctcrmined by other nunrcrical techniques. Wcrc. the cxistcncc of all these solutions is mentioned of’ cnhancinp the richness of’ the I-csptmsc occurring in the superharmonic
with the aim range ot’ the
systm.
3. ~‘ONCl,USIONS This study was conccrnecl with the nonlinear and chaotic response of R harmonically forced oscillator with quadratic and cubic nonlincaritics representative of ;I prohlcm in elastic structural dynamics. The results obtained show how the predicticm of irregular motions of a system (type, number, periodicity anti/or ‘chaoticity’) can bc made in two steps: (i) by obtaining stability boundaries of dominant approximate steady-state solutions-even of low order---in given ranges of values of the control parameters, e.g. around meaningful rcsonancc conditions; (ii) within unstable rcpions, hy identifying the solutions and looking for possihlc chaos through various dynamic measures on the basis of localized point-by-point and cell mapping computer simulations.
‘l‘hc
capabilily
situations
role
of
predict
satisfactorily
shown.
was
of computer
IW strcsscd. illilslrulc
An
It itllows: to
clucsliori
irregulnr
motions.
of
since
view,
unpredictable
I 7 _.
3. 4.
5. 0.
7.
x. Y.
LO.
1I. 12. 13. 14. 15. IO, 17 IX IY 20
In
sonic
(i)
kchnical
cases,
there
ilS
intcrcst with
to
the
iI is fcclt thal
roulcs
very
primarily
10 the expected
enc.
in the
dynamics;
of
I’rom
simple
‘l’hc
is again
to
kinds of strange
to initial
chaos
important
with
response
SystCill
sensitivity
Icrniination
is no(
to
in tliI’fcr_cnI
apl”“ximatiolls.
and diffcrcnl 01‘
strong
soluliow
occur
Cilll
of the sys~cm
IIlCiiSlIrCS
ZlCtUiil
this
itpprc)ximiitc chaos
higher-order
to chaos
rhcir
and
is concerned
respect
whcrc
clcscriplirm
rlU:iIltilittiVc
?lttlJCtOrS
(simple)
is ncccl for
differenI
rcliilhlc
of
regions
ii thorough
llli~lly
rcniains
I lowcvcr,
response
dclimilalc
10 identiiy of
arlitlysis
illld
for
CVi~lll~~lC
Chc cocxislcrice
open
numerical
simulations,
(ii)
iIllri~C(OrS;
stability bifurciitiofls
(iii)
10
conditions.
8 scqucncc a practical
01’
poinl
ass:~‘ssn~cn~ of
(he
Sohtonr & Fractais in Great Britain
Vol. 2, No. 3, pp.303-321,
1992 0
0960.0779/92$5.00 + .cPJ 1992 Pergamon Press Lrd
Prediction of Bifurcations and Chaos for an Asymmetric Oscillator F. BENEDETTINI, Dipartimento
di Ingegneria
delle Strutture,
Elastic
G. REGA and A. SALVATORI
delle Acque
e de1 Terreno,
Universitl
dell’Aquila,
L’Aquila,
Italy
Abstract-The nonlinear chaotic response of a harmonically forced elastic oscillator with quadratic and cubic nonlinearities is studied. The possibility of obtaining reliable prediction of bifurcations and chaos through combined use of stability analysis of low-order approximate solutions and accurate localized point-by-point and cell mapping computer simulations is examined. Some satisfactory results are obtained for the bifurcation predictive capability and the location of regions where chaos actually occurs.
1. BACKGROUND
AND MOTIVATION
Oscillators with asymmetric and symmetric nonlinearities have received increasing attention in recent years [l-7]. Besides the theoretical interest associated with the presence of both kinds of nonlinear terms, they have been specifically referred to in structural applications as representing the dynamic behaviour of elastic systems with initial curvature undergoing finite displacements. Systems governed by an equation of motion with quadratic and cubic nonlinearities in the form 4 + /Uj + c1q + c2q2 + c,q3 = PcosQt
(1)
have most often been studied. In Refs [4, 6, 71 a dimensionless oscillator of type (1) was considered as a one degree-of-freedom model to describe the planar finite dynamics of a prestressed suspended cable subjected to vertical harmonic forcing. Its coefficients (cl = 1, c2 = 35.953, c3 = 534.53) were given values representative of the properties of a cable used in overhead transmission lines and characterize it as a strongly nonlinear oscillator with a single equilibrium position. In particular, it was shown how, notwithstanding the model’s simplicity, the presence of even and odd nonlinear terms gives rise to a rich and complex set of responses [7]. A response chart in the forcing parameter space (Q, P) was built through point-by-point time integration of the motion equation with fixed initial conditions (Fig. 1). Various periodic solutions occur in the three main zones into which that space can be subdivided, namely the neighbourhoods of the l/2- and l/3-subharmonic resonance conditions (associated with the two kinds of nonlinearities, respectively), and that of the superharmonic resonance conditions. In each zone, with proper values of the forcing amplitude, the periodic solutions bifurcate to chaos following various routes and different kinds of strange attractors occur. The actual chaotic nature of those responses was tested using various qualitative and quantitative dynamic measures [7]. From the theoretical point of view, there is strong interest in obtaining reliable predictions of the conditions beyond which the periodic orbits are destroyed and persistent, steady-state strange attractors are established [8-lo]: unfortunately, such an objective is hard to obtain, so practically all developed criteria lead just to necessary conditions for 303
chaos. This is connected with the strong sensitivity of the response to control parameter variations and with the associated fractal nature of the basin boundaries not only in the initial condition space but also in the system parameter space. Thus, it is always possible to find subsets of parameter values leading to responses which look like irregular ones but cannot be directly classified as chaotic. However, such a circumstance is not a great drawback from the point of view of engineering system applications, for which the simple assessment of unpredictability of the response-or the onset of a dynamic regime different from the ordinary one-is often more important than its precise characterization as chaotic. With these problems in mind, it seemed worthwhile to examine the possibility of predicting the disappearance of the periodic motions obtained for system (1) and of locating limited regions in the system parameter space where irregular motions and chaos can occur, while concurrently analyzing the effect of variations of the i.c. Among the various predictive criteria-either based on the homoclinic orbit bifurcation [S-lo] or of heuristic type [9, 11,121 or associated with the stability analysis of regular nonlinear oscillations [l, 3,5,6, 13-16]-the latter were selected according to suggestions first made by Szemplinska-Stupnicka [l]. Indeed, they are of special interest in the present case because they allow the prediction of bifurcations which are possible precursors to chaos. and the establishment of relationships between the irregular and the dominant regular solutions occurring in various resonance zones. With these criteria, approximate steady-state solutions for the system are first obtained through one of various analytical methods and then the stability of the solutions is evaluated. One major procedural problem lies in the choice of the type and number of harmonic terms to be used for obtaining both accurate solutions and stability analysis. Approaches involving low [l, 3.5,6,13] or high-order [14-161 approximations can be distinguished, the latter are able, in principle, to furnish better results. However, two points are considered: (i) the main advantage of developing and applying criteria of this type is the possibility of recovering the huge patrimony of knowledge developed in the classical field of nonlinear oscillations for completing the understanding of new dynamic phenomena; (ii) interest is focused on identifying zones of response unpredictability rather than strictly chaotic motions. Within this ambit, it is worth examining how much accuracy can be obtained in the prediction of irregular responses by means of analytical solutions that are as simple as possible. Thus, low-order approximations will be considered later. Use of higher approximations is likely to be necessary in the study of more complex systems or phenomena, mostly as far as reliable stability evaluations are concerned [16], but these can sometimes be useless for engineering purposes due to roughness or uncertainty of the model considered. In the next section, the bifurcation predictive capability of approximate solutions properly chosen in the various regions of resonance is analyzed with the help of the results of accurate but localized computer simulations. 2. APPROXIMATE
PERIODIC SOLUTIONS, STABILITY ANALYSIS, BIFURCATIONS, OF ATTRACTION
BASINS
Equation (1) being asymmetric, its period jT solution q(t) does not satisfy the symmetry property, i.e. it is q( t + jT/2) # -q(t), j = 1, 2, 3, . . . . If one seeks a solution represented by finite Fourier series, an asymmetric solution generally contains a constant term and both odd and even harmonic components:
Prediction
of chaos for an asymmetric
305
oscillator
(2) can be determined through, e.g. the harmonic The coefficients of/j, ~i~j, in equation balance method. Examination of local stability of the steady-state period jT solution is made by adding to it a small disturbance, g(t) = q(t) + 6q(t). Then, inserting it into equation (1) and neglecting nonlinear terms in 6q(t), a linearized Hill’s type variational equation is obtained, 64 + ,LL& + g(q)dq
= 0,
(3)
in which the periodic coefficient g(q), evaluated through equation (2), has period jT. The Floquet theory [17] calls for solutions to this equation of the form (Floquet multiplier) of the monodromy 6q(t + jT) = A6q(t), where h. is an eigenvalue matrix C associated to a fundamental matrix solution Q, of equation (3), @(t + jT) = C Q(t). By assuming Q(O) = I, it results C = @(jT). The matrix Q(jT) can be computed by numerically integrating equation (3) for fixed values of P and S2 in the interval [0, jT] subject to two sets of i.c., (i) 6q, (0) = 1, 64, (0) = 0 and (ii) 6q, (0) = 0, 84, (0) = 1. The values of i\. determine the stability of the approximate solution q(t). IAl< 1, which means eigenvalues of C remaining inside the unit circle in the complex plane, is a necessary When varying a control parameter of the system, the and sufficient stability condition. manner in which d leaves the unit circle characterizes the local qualitative bifurcations + 1, one has occurring to the periodic orbit. When it leaves the unit circle through 6q(t + jT) = 6q( t). Thus, q(t) is a jT-periodic attractor and the system undergoes a jump bifurcation. On the other hand, when A leaves the unit circle through -1, it results 6q(t + 2jT) = as(t). Thus q(t) is a 2jT-periodic attractor and the system undergoes a period-doubling bifurcation.
2.1.
Neighbourhood
of l/2-subharmonic
resonance
The computer simulations of equation (1) with zero i.c. made to obtain Fig. 1, showed that the route to chaos from the right in the neighbourhood of the 1/2_subharmonic resonance is characterized by a dominant sequence of p.d. bifurcations [4,7]. This is clearly observed in Fig. 2, where the period chart of the response obtained at a fixed value of forcing amplitude (P = 0.04) is reported. To examine this phenomenon within the framework of stability analysis of approximate solutions, the underlying dominant period 1 response is described by the simplest asymmetric expansion, containing just the drift and fundamental harmonic: q(t) = a0 + alfcos(Q2f Substituting it into equation (1) and applying system of three algebraic nonlinear equations: 3c,(a2,alf
+ afJ4)
+ 8).
the harmonic
+ 2czaoalf + alf(l
(4) balance
+ 3a,a$)
leads
to a
- Q2) = Pcos6
alfpS2 = - Psin6 c,(2ai
method
(5)
+ c2(2ai + aff) + 2a. = 0.
By solving it using a numerical procedure, the unknown amplitudes and phases of the solution are obtained; in particular, the frequency-response curve can be built. Analysis of its local stability is made following the procedure outlined previously. The eigenvalues of C have been repeatedly computed as Q is varied at fixed values of P. The frequency-reszones of instability (thin lines) are ponse curve for the amplitude aIf and the relevant plotted in Fig. 3(a) (P = 0.04). Besides the classical unstable hysteresis branch, a further
306
F. BENEDETTINI et al, 1.0
0.0
1’
0.4
0.4
0.3
0.3
0.2
0.2
0.1
0.1
0.0
---1
1.0
0.0
I-3
Period
Period
1 m
6xN
Fig. 2. Period
chart
1.70
I
3.0
Period
m__
Fig. 1. Periodic
1.66
I
2.0
1.74
2 i
Period
and chaotic
I .7a
in the neighbourhood
response
1.82
3
in control
1.86
of l/2-subharmonic
R
Period
m
Transition
0.0
I
4.0 6
Chaos space (P, n).
1.90
resonance
1.94
1.98
(P = 0.04, zero i.c.).
Prediction
of chaos for an asymmetric
oscillator
307
large zone of instability is attained through p.d. bifurcations on accordingly, the establishment of a period 2 solution is expected. Thus one can consider the simplest expansion of period 2, namely q(t) and obtain amplitudes
= a0 + a,cos(Qr
+ 79) + a,,,cos(Q/2~t nonlinear method,
,/+os(6 ($3% + $4Q2
+ 2C,U” + &z:
- a,(@ (ic3u0 (3c,u,u,
+ ic2)a&sin(t9
+ c2ul)ul/2COS(19
-
lower
branch:
+ q),
the following system of five algebraic and phases through the harmonic balance - 299 + a,(3c3a;
the
equations
(6) in the
unknown
+ ;C&)
1) = P cos 6
- 2cp) - a,,~&? = Psinff
- 2q) + 2c*a~~~,,/ +
c3
i
3u; + ; u: + ; a:,2 a,/2 31
!a2 + 1-qu:/2=0
1
i (3c3uoul ~c3u&c0s(6
+ c~u~)u~,~ sin (79 - 2y) -
+LZ~/~@ = 0
- 2l7I) + c*[ui + ;
and solve it numerically. At the considered value of forcing amplitude, the range of existence of this solution coincides precisely with the range of second-order instability of the period 1 solution: the frequency-response curve for the amplitude ali of the subharmonic component is plotted in Fig. 3(a), as well. By repeating the stability analysis, the lower branch is seen to be practically always stable, while the upper branch shows a large zone of instability attained again through p.d. bifurcations: the occurrence of a period 4 solution is thus expected.
(4
(b)
cl.0
0.5
1.0
Forcing Fig. 3. Frequency-response
15
2.0
2.5
3.0
3.5
4.0
frequency curves: approximate solution (thick, stable and thin, unstable) lation results with own i.c. (black, also zero i.c.). (P = 0.04).
and computer
simu-
308
F. BEN~D~TINI et al
Use of higher subharmonic approximations would make likely the possibility of further p.d. bifurcations evident. However, these occur in very small ranges and the desire is to limit the analysis of the bifurcation predictive capability to low-order approximations. So, the predictions obtained are compared directly with the results of computer simulations. The amplitudes of the harmonic and 1/2subharmonic components obtained by calculating the FFT of the simulation results are compared in Fig. 3(b) with the corresponding harmonic balance curves. Different from the stability analysis, numerical simulations require specification of the i.c. The values q(O), Q(O), obtained by inserting the harmonic balance amplitudes and phases in equations (4) and (6), are considered here (white and black dots in Fig. 3(b)) compared to the zero values already ascertained for obtaining the chart in Fig. 1 (black dots in Fig. 3(b)). It is observed that by using the conditions defined by equation (4), the upper a,, branch is almost completely matched, while its lower branch is matched only where it is stable according to the approximate analysis. This is not shown in Fig. 3(b). Where the analysis indicates occurrence of p.d. bifurcations, the conditions allow the curves a,‘?, a, of the period 2 solution to be followed. When using i.c. defined by equation (6) and decreasing the frequency, the lower u,/? branch is matched from the frequency value (Q G 2.26) where the period 2 solution arises up to a point close to the vertical tangent. This is consistent with the branch’s observed stability and, thus, agreement with the theoretical curve is very good. The response then jumps to the upper ali branch. The upper branch is followed from its peak amplitude up to just the point where the stability analysis shows the occurrence of a p.d. bifurcation. Left of this point, a l/4-subharmonic component appears in the numerical solution and then successive p.d. bifurcations to period 8, 16, . . are observed. As the bifurcation proceeds towards chaos, the amplitudes of 1/2_subharmonic component decrease slightly with respect to the theoretical values, due to the presence of further subharmonic components. The zone of chaos extends approximately between Q = 1.82 and Q = 1.70. Below this value, even with the i.c. defined by equation (6), the response jumps to upper period 1 and the relevant ~7,~ branch is followed up to a frequency value (Q z 1.52) where the response jumps down to period 2 and the corresponding theoretical curves-which are now stable again-can be followed numerically. The analysis confirms substantially the points made in [I], where the same problem was studied for a system with asymmetric nonlinearity similar to equation (1) by considering different values of the system parameters. Slight differences are represented by the occurrence of the chaos zone between the two points of vertical tangencies and by the possibility of moving along the al:? branch left of that zone: the latter suggests that the irregular motions found are actually associated with the specific second-order instability of the period 2 solution. When decreasing the frequency with zero i.c., the lower a1;1 branch is followed up to a value of frequency (Q = 1.96) below which the process of bifurcation to chaos develops and solutions with higher period (e.g. period 4 at B = 1.92) are established exhibiting the 1/2subharmonic component with upper amplitude. In substance, the jump to the upper branch has shifted to a slightly higher value of frequency with respect to the previous case. These local shifts heuristically suggest the possibility of a fractal-like nature of the boundaries between different periodic (and/or chaotic) regions in the forcing parameter space depending on the i.c. The capability of stability numerical analysis of the approximate solutions to accurately predict bifurcations and delimitate regions where chaos can occur is illustrated more limit curves of generally in Fig. 4 in the forcing parameter space. Here, the stability approximate period 1 (equation 4) and period 2 (equation 6) solutions and the locus of
Prediction
of chaos for an asymmetric
309
oscillator
points of vertical tangency to the period 2 frequency-response curve are superimposed on the regions of response obtained through numerical simulations around the l@-subharmanic resonance. Of course, for a correct comparison of numerical and theoretical results, one must remember that the former refer to zero i.c. These belong in various zones: to the domain of attraction of the lower period 1 solution for frequencies higher than the corresponding right stability limit curve (zone 1 in Fig. 4); left of this curve, to the domain of lower period 2, that is always stable (zone 2); where lower period 2 does not exist, to the domain of upper one (zone 3); finally, to the domain of upper period 1 which is always stable (zone 4). Embedded in these zones and even nearly bounded on the right by the limit curves of first or second order instability of period 2 solution, a large zone of irregular (highly periodic) and chaotic motions occurs. To clarify all features of the system dynamic behaviour in the considered range of forcing parameters, one has to complement the regions of periodic and chaotic response of Fig. 4 with systematic construction of basins of attraction in the i.c. space. In particular, the extension of the domain of chaotic solutions must be checked. Shown in Fig. 5(a-f) are the basins of the different solutions obtained with a cell mapping technique [18] for various values of decreasing frequency at the same forcing amplitude as in Fig. 3. To the right of the zone of existence of the period 2 solution, the basins of the upper (light blue) and lower (blue) period 1 spiral into each other. At the value Q = 2.26 (a), corresponding to the peak amplitude of the ali curve, very thin tongues (yellow) relevant to the upper period 2 solution appear within the basin of lower period 1. Just to the left of that value, consistent with the p.d. bifurcation prediction of this latter solution, the blue basin corresponds to lower period 2, though preserving the same knob shape. Moving towards lower frequency values, the two basins of lower and upper period 2 solutions evolve differently: the former one vanishes while the upper one increases, though
l .
.
l
Stability
limit
of
period
1
solutron
o
Stability
limit
of
period
2
solution
-
Limit
of
vertical
tangency
in
sub.
l/2
i 0.1
i
Period
El
, Ij ‘$
----I I!zrI
Transition
Penod
1.3 Fig. 4. Bifurcation
2.0 predictive
capability Comparison
1
to
3
Period
i
2
chaos
Chaos
2.6 of stability numerical analysis of period 1 and 2 approximate with computer simulation results (zero i.c.).
solutions.
310
F. BENEDETTINI et al
corresponding to solutions with higher and higher periodicity which eventually end in chaos (b-d). The strong fractal nature of the relevant basin boundaries is worth observing (b, c). At the same time, the size of the knob diminishes and the basin of attraction of upper period 1 solution increases notably. The ‘sudden’ change of the response from chaotic to upper period 1 occurring at the left boundary of the chaos zone in Figs 3(b) and 4. corresponds to strong mixing of the relevant basins in the i.c. space, that however is observable only through use of a much narrower window (e). To the left of that zone, any i.c. leads to the upper period 1 solution up to the frequency value (Q z 1.52) where a new spiraling tongue, corresponding to solutions with higher periodicity, appears in the i.c. space and the period 2 solution is again established through a reverse p.d. bifurcation (f). As the frequency is further decreased, this tongue becomes thinner and thinner and when the period 2 solution disappears, it corresponds to the lower period 1 solution still occurring in a short frequency range (see Fig. 3(b)). The important role played by the i.c. clearly emerges. Moreover, it is worth noting how the observed instability of the dominant approximate periodic solution-though a necessary condition for chaos-also suggests the most likely i.c. to consider for actually obtaining it. In the present case, they were rather successful.
2.2.
Neighbourhood
of l/3-subharmonic
resonance
subharmonic and chaotic responses According to the results of computer simulations, arise in this frequency range following more mingled routes with respect to those occurring in the neighbourhood of l/2-subharmonic resonance. This can be observed in the period arise at certain frequency chart of Fig. 6, again obtained with zero i.c. Period 3 solutions values and intertwine with period 6 solutions into which they bifurcate themselves. The latter give rise to several responses with period multiple of 6 increasing irregularly up to very high values. Accordingly, the chaotic strange attractor exhibits six independent bundles orginating from the period 6 solution around which the motion fluctuates [7]. Besides the fundamental solution equation (4), approximate steady-state solutions are sought by considering the simplest expansion of period 3, namely q(t)
= a,, + a,cos(Qt
+ 79) + a,/,cos(Q/3-
f + v).
_.~~_
1
Ch.
4
Fig. 6. Period
chart
in the neighbourhood
of l/3-subharmonic
resonance
(P = 0.40, zero i.c.)
(8)
Prediction
of chaos for an asymmetric
oscillator
Fig. 5. Basins of attraction in the neighbourhood of l/2-subharmonic resonance. P = 0.04; P = 2.26 (a), 2.00 (b), 1.92 (c), 1.76 (d), 1.69736 (e), 1.50 (f). Range for q(O), Q(O): -0.2 + 0.2 (-0.02 + 0.02 in case (e)).
312
Fig. 8. Basins
F. BENEDE~INI
etal
of attraction in the neighbourhood of l/3-subharmonic resonance. P = 0.40; B = 3.35 (a), 3.42 (b), -0.3 + 0.3; for Q(0): -0.4 + 0.4. 3.58 (c), 3.82 (d), 4.13 (e), 5.00 (f). Range for q(0):
Prediction
of chaos
for an asymmetric
oscillator
313
0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2 0.08
0.08
0.07
0.07
0.06 g g 0.05
;..
&.. . 8004 ! _4 zi g ................_./_..._....._ .. ..
I . . . . . . .
...)
eo.03
‘_ j
.
.......
0.06
\ ......
0.05
(
......
0.04
,
......
0.03
2
0.02
.
.
. .
: ......
.
0.01
0.02
0.01
0.00 0.00 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.3 0.3 1.0 1.1 1.2
"A,
EA,
Fig. 11. Frequency-response
Xi-A,OA, A,*A6.A,~A,AA,
Forcing frequency
curves: approximate solution (thick, stable and simulation results with own i.c. (P = 0.10).
thin,
unstable)
and
computer
314
F. BENEDE~INI
et al.
Pig. 13. Basins of attraction in the neighbourhood of order 2 superharmonic (b), 0.89 (c), 0.90 (d), 0.91 (e), 0.93 (f), 0.954 (g), 0.95526 (h). Tl green,
peak. P = 0.10; Q = 0.85 (a), 0.87 T2 yellow, T3 light green, T4 blue.
Prediction
of chaos for an asymmetric
By applying the harmonic balance method, in the amplitudes and phases is obtained 14 c&c0s(6
- 37#) + a,(3c,a;
2C$z~ + +ju:
+
i c,aT,3sin(7? 3 4
the system
315
oscillator
of five algebraic + +&j)
nonlinear
- a@
equations
- 1) = Pcos6
- 3rjJ) - alpQ = Psin6
3 - 3r/~) + 2c2a0al/3 + cx 3ai + 5 a: + 43 ali3 2)al,~+(l-~)&?l,3=0 (
c,alf&cos(6
: c3ala:j3sin c2[ai + ;(a:
(6 - 3q)
+ f a,ppQ
= 0
(9)
+ ~75~)] + c3[ai + tao(a: + ~z:,~)] + a0 = 0,
and then solved numerically. The frequency-response curves for the period 1 solution and for the amplitude ali3 of the period 3 solution are plotted in Fig. 7(a); a value of forcing amplitude (P = 0.40) is considered for which the simulations showed occurrence of chaos (Fig. 6). The stability analyses show that both the upper and the lower branch of period 1 are always stable in the range of interest, while the upper branch of the subharmonic component exhibits a large zone of instability (3.65 < Q < 4.98) attained through p.d. bifurcations: at the branching points, the occurrence of a period 6 solution is thus expected. Again, the predictions are compared directly with the results of computer simulations obtained with both zero i.c. (black dots) and the conditions corresponding to the harmonic balance approximations (white and black dots) in Fig. 7(b). With those defined by equation (4), the two stable branches of aIf are completely matched except in a very small part of the lower one near the vertical tangency to the aI/3 curve. With the conditions defined by equation (S), period 3 responses are obtained with very high accuracy in the whole range of theoretical existence of l/3-subharmonic solution, except in a narrow zone (3.60 < Q < 3.74, the same found when using zero i.c., see Fig. 6) where period 6 solutions appear and, correspondingly, the amplitude a,/3 decreases slightly with respect to the theoretical value. The responses obtained with zero i.c., already shown in Fig. 6, are also
(a)
(b)
0.0
2.0
Forcing Fig. 7. Frequency-response
40
6.0
frequency
a.0
10.0
IL
3.00
~~~~~
~
4.00
Forcing
5.00
6 00
frrquency
curves: approximate solution (thick, stable and thin, unstable) lation results with own i.c. (black, also zero i.c.).(P = 0.40).
and computer
simu-
F.BENLD~TTINI et al
316
indicated from the black dots and the arrow in Fig. 7(b), where the frequency range in which chaos occurs is also shown. The bifurcation from period 3 to period 6 response evident from previous results is consistent with the bifurcation prediction of the stability analysis of the period 3 approximation (8). However, the notably different kind of response (period 3 or chaotic) obtained with different i.c. raises some doubt about the overall reliability of the analysis and suggests a strong influence of the i.c. on the results. Better understanding of the complex set of responses which occur can be achieved by examining the evolution of the relevant basins of attraction with increasing frequency (Fig. 8(a-f)). To the left of the zone under investigation, only a small basin of lower period 1 (violet) occurs within upper period 1 (red) (a). Where the stable period 3 solution arises, thin tongues (yellow) enter the former basin (b), and then quickly expand over practically the whole basin. Thereafter, consistent with the prediction of p.d. bifurcation of the relevant stability analysis, this basin corresponds to solutions with increasing period (6, 12, . .) ending in chaos (c, d). Two points are worth noticing. First, the basin of attraction of lower period 1 is reduced up to the knob core, but it never vanishes: this is consistent with the overall stability of the solution shown by the approximate analysis. Second, the complex evolution to chaos of the response within the knob is accompanied by the onset of new tongues of period 3 solution (pink) surrounding the core of lower period 1 basin (d, e). Thus, the period 3 response obtained for its own i.c. within the chaos zone and charted in Fig. 7(b), is associated only with the existence of these tongues, to which those conditions belong. However, though taking into account the strong mixing of the various basins occurring within the knob, most of the relevant i.c. lead to responses with period other than 3 (lower periods 1, 2, 4 and 6 in Fig. 8(e)), thus confirming in general terms the instability prediction of that solution obtained with the approximate analysis. Such a circumstance is understood more clearly through comparison of the detailed picture of the response regions obtained with computer simulations in the zone of main interest around the l/3-subharmonic resonance, and the regions of stable and unstable period 3 solutions given by the approximate analysis (Fig. 9). A shift of the approximate unstable region towards higher frequencies with respect to the numerical one is observed. But the overall picture is very similar in the two cases, with the wedge of stable solutions entering the upper part of the figure. Moreover, the shift has limited meaning since the numerical results correspond only to zero i.c. On the other hand, the basins of attraction obtained at different forcing amplitudes show the period 3 solution tends to occupy the whole knob in the upper part of Fig. 9, consistent with its approximate stability. Instead. it is greatly reduced by the occurrence of various further solutions-with marked fractal character of the knob itself-in the lower part of Fig. 9, where the period 3 solution practically vanishes (e.g. it does not exist at P = 0.2, Q = 3.25). Finally, when considering the even greater role played by the i.c. in this parameter zone, the capability of the simplest approximate solution to delimitate regions where complex and chaotic motions can occur is satisfactorily established. Coming back to Fig. 8, it is observed that with the further increase of the frequency (f), the basin of period 3 solution (yellow), now stable again, increases notably and remains the only one in the knob-together with the lower 1 basin-while the knob itself tends to occupy the whole window.
2.3.
Neighbourhood
of superharmonic
resonances
the response picture is quite rich and complex with To the left of primary resonance, many periodic and chaotic zones not so well defined, but following a certain regular pattern
Prediction
of chaos for an asymmetric
317
oscillator
. .
0.46
.
1 .
.
limit 3
Limit of vertical tangency in sub
-
.
Stability of period solution
. .
0.36
3 .3 Fig. 9. Regions
3.5
3.7
of stable
3.9
4.3
4.1
and unstable
4.5
4.9
4.7
5.1
o
Period
1
o
Period
3
[1
Period
1,3
D
Period
1,3,6
u
Period
6
o
Period
6*N
-
Chaos
5.3
period 3 solutions with approximate simulation results (zero i.c.).
analysis.
Comparison
with computer
of inclined and nearly parallel stripes embedded within the dominant period 1 solution and extending over the entire range of forcing amplitude (Fig. 1). Several transition zones occur in which the response periodicity is notably sensitive to small variations of forcing parameters (see also the period chart of Fig. 10, still obtained with zero i.c.). It is thus easy to predict high sensitivity to variations of i.c., too. An approximate steady-state solution is sought by considering a period 1 solution improved by the presence of higher harmonics, namely the order 2 and 3 superharmonic terms: + p) + a,cos(3Qt
q(r) = a(J + UlCOS(Qf + ?Y)+ azcos(251t
__
Chaos
_._
+
I@).
_--I-I_
-
72 a
0
‘f: g
20
Qi z k !.G p:
aL 4 1
-_A
0.2 Fig. 10. Period
.:--y.:---: .. . .
chart
0.3
.. .
__--
_. ___
0.4
in the neighbourhood
.0.5
‘t_
0.6
of superharmonic
_
-
R
-
0.7 resonances
---~ 0.8
0.9
(P = 0.10, zero i.c.).
(10)
318
F. BENEDETT.INIet al
With the usual technique, 3c+ 4
[afcos(38
the system
- I/J) + a:cos(6
of seven algebraic
nonlinear
equations
- 2p + q)] + az(cz + 3c,a,)[a,cos(26
reads:
- cp)
3c3 a: + fp - q)] + a, ao(2c2 + 3cjao) + y y + ai + ai i [ I
+ a3cos(6
1 43c3a3
= Pcos-9
+ (1 - n*)
[afsin(38
-
+ a,sin(@ al(c2 + 3c,a,)
q) + azsin (6 - 2~ + q)] + a2(c2 + 3c,ao)[aI sin (2tY - q) q)] - ,@a,
+ q -
- q) + a3 cos(B + cp - tp> + 3c,a,
$cos(211 I
1 al(c2 +3c3ad [ + ai + 2
?
+ 2
7-
ai + 4
Tsin
-
3c3ala2a3
2
+ a,(1 - 4Q2 + 2c2a,)
(21) -
__
2q1-t
l/Q)
= 0
q) + a3 sin (8 + q -
1
(11)
ye)
sin(8-2fp++)-2yQa2=0
3
7cos(31’)
ala3 1 I2cos(6-
= Psint9
- v) + ala2 (c2 + 3c,a,,)cos(6 i
+ i 22 3c,a3 ai + :
+ :
3c3a7 + y - I/I) + -ycos(6
- 2y, + q)
ai + 4 i + a,(1 - 9R2 + 2c,a,)
1
= 0
3 ysin(3B
3c3ala2
-
ijj) + ala2 [ (c2 + 3c,a0)sin(6
+ 3pQa3 = 0 1 ~cos(213-cj7)+~cos(r’l+~-~) + a0(c2aC, + c,ai
+ q -
1
v) - Tsm(8 3c3a2
i +(cZ+3c3a,,)~+$+$
2
- 295 + I./I)I 2
2
1
+ 1) = 0.
The frequency-response curves obtained for a,, a2, a3 at the same value of forcing amplitude as in Fig. 10 are plotted in Fig. 11. The considerable modification of curve a, of the fundamental harmonic associated with the two superharmonic resonances Q = l/2 and Q = l/3 is clearly observable. Correspondingly, the superharmonic amplitudes a2 and a3 increase strongly [19,20]. The construction of the curves through numerical solution of system (11) was not forced upward to obtain the unstable branches of the loops which are likely to occur between the points of vertical tangency for both the harmonic and the superharmonic components, the former downward and the latter upward [19], but they are strongly suggested by the trends of the curves themselves. The stable parts of the approximate solution are denoted as usual with thick lines in Fig. 11, but these are limited only to the a, curve for the sake of clarity. In the same figure the results obtained from computer simulations, made on a very fine frequency mesh by
319
Prediction of chaos for an asymmetric oscillator
considering i.c. corresponding to harmonic balance approximation (lo), are also plotted. The overall agreement is satisfactory. Discrepancies occur: (i) as far as the solution is concerned, mostly in the lowest frequency range, where higher superharmonics (namely, 4 and 5)-obtained from computer simulations and also plotted in Fig. 11-play an important role; (ii) as far as stability is concerned, with respect to the localization of some unstable zones within the frequency range. Consistent with the findings from Fig. 10, the picture of motions ensuing from the basins of attraction is very complex. Thus it is not possible to follow their evolution, but only to focus on some local results. Three main frequency zones can be distinguished in Fig. 11, two located in the neighbourhoods of the order 3 and 2 superharmonic peaks, the third one just between them. As to the first zone, it can be seen that to the left of the peak (Q G 0.52), the i.c. corresponding to both upper and lower approximate period 1 solution belong to the upper 1 basin: this is consistent with the stability and instability predictions obtained respectively for the two solutions (Fig. 11). By slightly increasing the frequency, a chaotic basin-to 1 belong-spirals within the upper period 1 basin which the i.c. of lower period (Fig. 12(a)), that decreases (b) till the solution itself vanishes (vertical tangency at Q G 0.5805 in Fig. 11). At the same time, the solution in the growing basin-which eventually remains alone-evolves quickly and irregularly through many different degrees of periodicity. In the second zone, corresponding to the intermediate branch of fundamental harmonic in Fig. 11, a rather variable response picture occurs in the i.c. space, with one or more basins of low periodicity. Between about Q = 0.7 and 0.8, corresponding to a decrease of basin fills in the overall i.c. window. the simulation amplitude al, a chaotic or subharmonic It is worth observing how the approximate analysis reveals a zone of unstable period 1 solutions, but this is shifted towards higher frequency values with respect to the simulation results. Accounting for higher superharmonics in the approximate solution would be likely to be necessary for better stability prediction [16]. The third zone corresponds to the order 2 superharmonic peak and is quite interesting. The basins of attraction with the increasing of frequency (Fig. 13) evolve from a single
(b)
Fig. 12. Basins
of attraction
in the neighbourhood (b). Range
of order 3 superharmonic for q(O), 4(O): -0.2 + 0.2.
peak.
P = 0.10;
Q = 0.56 (a), 0.58
period I basin (?‘I) to four coexisting hasins (;I-d)--- all of period I --appearing in (urn (‘1’2. ‘1’3, ‘1’4). For SZ ~~0.90. the nunihcr of period I solutions diminishes down to one’ (c-h) hut with ;I vanishing scqucncc different from the appcarancc scqucncc. Indeed, the surviving basin corresponds to ;I motion (‘1‘3) other than the original enc. ‘I’hcsc two solutions must Iw rclwtcd to the left and right hranchcs of the I’untla~nc~~t~~l in f;ig. I I. A clear undcrstanclil~g of what happens is xhicvctl
h;~nnonic through
(I, around the peak a comparison ol’ I:igs
I3 and 14. The latter shows an cnlargcmcnt of the 0, curve around the peak, ohtainctl from computer simulations with proper i.c. 13esidcs the two main period I branches--‘I’1 from the left and T3 to the right- which overlap and arc likely to co~~ncct with WC‘ another below via an unstahlc loop 1201, two more period I brnnchcs (‘I’2 and ‘1‘4) arise in ;I small frequency range and exhibit nearly the same (7, amplitutlcs as the previous ones. ‘Thcsc two branches overlap in a short range as well. where four stahlr solutions exist, and at-c’likely to connect with one another via unstable branches that could he dctcrmined by other nunrcrical techniques. Wcrc. the cxistcncc of all these solutions is mentioned of’ cnhancinp the richness of’ the I-csptmsc occurring in the superharmonic
with the aim range ot’ the
systm.
3. ~‘ONCl,USIONS This study was conccrnecl with the nonlinear and chaotic response of R harmonically forced oscillator with quadratic and cubic nonlincaritics representative of ;I prohlcm in elastic structural dynamics. The results obtained show how the predicticm of irregular motions of a system (type, number, periodicity anti/or ‘chaoticity’) can bc made in two steps: (i) by obtaining stability boundaries of dominant approximate steady-state solutions-even of low order---in given ranges of values of the control parameters, e.g. around meaningful rcsonancc conditions; (ii) within unstable rcpions, hy identifying the solutions and looking for possihlc chaos through various dynamic measures on the basis of localized point-by-point and cell mapping computer simulations.
‘l‘hc
capabilily
situations
role
of
predict
satisfactorily
shown.
was
of computer
IW strcsscd. illilslrulc
An
It itllows: to
clucsliori
irregulnr
motions.
of
since
view,
unpredictable
I 7 _.
3. 4.
5. 0.
7.
x. Y.
LO.
1I. 12. 13. 14. 15. IO, 17 IX IY 20
In
sonic
(i)
kchnical
cases,
there
ilS
intcrcst with
to
the
iI is fcclt thal
roulcs
very
primarily
10 the expected
enc.
in the
dynamics;
of
I’rom
simple
‘l’hc
is again
to
kinds of strange
to initial
chaos
important
with
response
SystCill
sensitivity
Icrniination
is no(
to
in tliI’fcr_cnI
apl”“ximatiolls.
and diffcrcnl 01‘
strong
soluliow
occur
Cilll
of the sys~cm
IIlCiiSlIrCS
ZlCtUiil
this
itpprc)ximiitc chaos
higher-order
to chaos
rhcir
and
is concerned
respect
whcrc
clcscriplirm
rlU:iIltilittiVc
?lttlJCtOrS
(simple)
is ncccl for
differenI
rcliilhlc
of
regions
ii thorough
llli~lly
rcniains
I lowcvcr,
response
dclimilalc
10 identiiy of
arlitlysis
illld
for
CVi~lll~~lC
Chc cocxislcrice
open
numerical
simulations,
(ii)
iIllri~C(OrS;
stability bifurciitiofls
(iii)
10
conditions.
8 scqucncc a practical
01’
poinl
ass:~‘ssn~cn~ of
(he