Non-linear dynamics of multiple friction oscillators

Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

www.elsevier.com/locate/cma

Non-linear dynamics of multiple friction oscillators Ugo Galvanetto 1 Dipartimento di Costruzioni e Trasporti, Universit a di Padova, Via Marzolo 9, 35131 Padova, Italy Received 28 July 1998

Abstract This paper deals with the dynamics of multiple friction oscillators. In some particular cases, when the driving velocity is small, the dynamics of the systems can be described by low-dimensional discrete maps which allow all existing attractors and their basins of attraction to be detected. Moreover the maps allow the relevant Lyapunov exponents to be computed and the classical path following techniques to be applied. The complex dynamics of the systems is clearly illustrated. Ó 1999 Elsevier Science S.A. All rights reserved.

1. Introduction This paper presents some results concerning the non-linear dynamics of stick±slip systems composed of a chain of blocks moving on a rough surface. The dynamics of this kind of system is non-smooth since a generic motion is composed of a sequence of stick and slip phases in which the oscillator sticks on the rough surface with no relative motion between them or slips with respect to the same surface. Moreover, in stick± slip systems, some state variables are blocked during the stick mode; thus, they belong to the class of nonsmooth systems with a variable structure, where the degrees of freedom change in time. The dynamics of these systems is, therefore, governed by di€erent equations according to whether the stick or slip modes occur and its numerical integration requires some caution as explained in Ref. [1]. In the past, in the ®eld of applied sciences, stick±slip systems were investigated in seismology in particular because fault dynamics was often modelled by means of chains of blocks undergoing stick±slip motions on moving surfaces [2±5]. More recently, also in the ®eld of mechanical engineering renewed attention has been drawn to such systems by Popp and co-workers [6±8]. Finally the rich and complex dynamics of these models has been investigated by some authors in the spirit of non-linear dynamics and chaos [9±11]. This paper belongs to the latter group and aims at giving a more general description of the stick±slip dynamics of chains of blocks. In Section 2 the single stick±slip oscillator is presented to introduce the basic ideas which are then developed in Section 3 where the multiple friction oscillator is de®ned. In Section 4 a fundamental low-dimensional map is described and some examples for two- and three-block chains are also explained. Finally the results are shown in Section 5 and some conclusions are also given. 2. The single degree of freedom mechanical system The simplest example of a dry friction oscillator is shown in Fig. 1: a mass m, supported by a moving belt, is connected to a ®xed support by a linear elastic spring of stiffness k. The contact surface between

1

Tel.: +39-49-827-5585; fax: +39-49-827-5604; e-mail: [email protected]

0045-7825/99/$ ± see front matter Ó 1999 Elsevier Science S.A. All rights reserved. PII: S 0 0 4 5 - 7 8 2 5 ( 9 9 ) 0 0 0 2 1 - 3

292

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

Fig. 1. Single block dry friction oscillator.

mass and belt is rough so that a friction force is exerted by the belt on the mass. Energy is continuously introduced into the system by the continuous motion of the belt which moves at a constant speed vdr (called driving velocity) while the friction force may either introduce or dissipate energy according to whether the friction force has the same sign as the velocity of the oscillator or the opposite. This system is the prototype of mechanical systems undergoing self-sustained oscillations and is presented in classical textbooks [12±14]. The stick motion is uniform since the oscillator moves with the belt at the same velocity, this motion can be described by the law of the motion: x…t† ˆ vdr t ‡ x0 ;

…1†

where x, the displacement of the oscillator, is measured from the position in which the spring assumes its natural length, x0 is the initial displacement and t the time variable. As the displacement x of the block grows the elastic force increases up to the point in which it is equal to the maximum static friction force Fs . Therefore the position of impending slip is given by x ˆ Fs =k:

…2†

The slip motion is governed by an equation of the following type: ma ‡ kx ˆ ld …v ÿ vdr †Fs ;

…3†

where v ˆ dx=dt is the velocity, a ˆ d2 x=dt2 the acceleration and ld …v ÿ vdr † the dynamic friction coecient which is equal to 1 at impending slip when the relative velocity between oscillator and belt, vrel ˆ v ÿ vdr , is zero. Eq. (3) depends on several parameters which make its study complex, therefore an adimensionalisation procedure can be applied in order to reduce the number of parameters. If the following dimensionless variables are introduced in Eq. (3) r k k t; X ˆ x; …4† sˆ m Fs where s is the non-dimensionalised time and X the non-dimensionalised displacement, Eq. (3) becomes: X ‡ X ˆ Md …X_ ÿ Vdr †;

…5† p where the upper dots indicate derivation with respect to dimensionless time, Vdr ˆ vdr km=Fs is the dimensionless driving velocity and Md …X_ ÿ Vdr † the dimensionless dynamic friction force, which is a function _ _ of the value pof  the relative velocity Xi ÿ Vdr , equal to 1 at impending slip when X ÿ Vdr ˆ 0. Note that d=dt ˆ k=m…d=ds†. In the dimensionless form the impending slip condition (given by Eq. (2) in dimensional form) reads X ˆ 1. It is worth observing that the block may stick on the belt in the interval ÿ1 < X < 1, but out of this interval stick motion is not possible, such an interval is the locus of possible stick phases.

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

293

The dimensionless friction characteristic is assumed to be given by the following equation: Md …Y ÿ Vdr † ˆ

1ÿa 2 ‡ a ‡ b…Y ÿ Vdr † ; 1 ‡ cY ÿ Vdr

0 6 a < 1; b P 0; c > 0

…6†

which is illustrated in Fig. 2 for di€erent combinations of the relevant parameters. c is positive because the dynamic friction force has to be decreasing for small values of the relative velocity, b is positive because the dynamic friction force is usually assumed to be increasing for large values of the relative velocity, a varies between 0 and 1 to avoid unrealistic changes of sign in the friction force [6,15]. It is possible to demonstrate [9] that the only possible ®xed point of this mechanical system, …X ˆ Md …ÿVdr †; X_ ˆ 0†, is unstable and that the only possible attractor for the one-block system is a limit cycle [16], but the peculiar nature of stick±slip systems requires some additional considerations about the di€erent types of periodic motions. Having ®xed the frictional properties let us keep the driving velocity as the only variable parameter: (a) for a small driving velocity the limit cycle is shown in Fig. 3. In this case the velocity of the block is always less or equal to the velocity of the belt, for most of the time the block rides on the belt with no motion relative to it and short slipping phases separate relatively long stick phases;

Fig. 2. Friction characteristics for di€erent values of the relevant parameters: (a) a ˆ 0, b ˆ 0, c > 0, (b) a > 0, b ˆ 0, c > 0, (c) a > 0, b > 0, c > 0.

294

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

Fig. 3. Attractor for low driving velocity, Vdr ˆ 0.16: (a) phase space, (b) time history, (c) velocity in time, (d) friction force in time.

(b) for a medium driving velocity the limit cycle is shown in Fig. 4. In this case the velocity of the block may be less, greater or equal to the velocity of the belt, the stick phase is much shorter than in the previous case; (c) for a high driving velocity the limit cycle is shown in Fig. 5. In this case there is no stick phase and the block keeps slipping assuming velocities which are higher or lower than the velocity of the belt and only instantaneously equal to it. These ®gures also show that the system is subject to discontinuous friction forces. Our study is limited to the stick±slip cases (a) and (b) which are stick±slip motions. Fig. 6 shows the periodic attractor of the dynamic system described by the parameter values: a ˆ 0, b ˆ 0, c ˆ 3, Vdr ˆ 0.11. The segment ÿ1 < X < 1 on the line X_ ˆ Vdr is the geometrical locus of the stick motions: if velocity and position of the block belong to that segment the block sticks on the belt up to the point X ˆ 1 where a

Fig. 4. Attractor for medium driving velocity, Vdr ˆ 1.6: (a) phase space, (b) time history, (c) velocity in time, (d) friction force in time.

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

295

Fig. 5. Attractor for high driving velocity, Vdr ˆ 16.0: (a) phase space, (b) time history, (c) velocity in time, (d) friction force in time.

Fig. 6. Any initial condition of the phase space gives rise to a motion which exactly converges onto the attractor in a ®nite time.

296

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

permanent slipping phase starts. All the points on that segment give rise to motions which eventually coincide with the periodic attractor. The transient trajectories in the regions 1 of the ®gure arrive at a point of the previously mentioned segment and therefore are exactly attracted on the attractor in a ®nite time. The transient trajectories in the regions 2 of the ®gure end up in a point of the regions 1 and therefore they also are exactly attracted on the attractor in a ®nite time. The same argument applies to any point of the phase plane, consequently all initial conditions of the phase plane at a ®nite distance from the attractor converge on it in a ®nite time. The numerical investigations reveal that this is also true for the region of the phase space internal to the periodic cycle (with the exception of an in®nitesimal neighbourhood of the ®xed point …X ˆ Md …ÿVdr †; X_ ˆ 0†. Therefore, loosely speaking, we can assume that the largest Lyapunov exponent of the periodic stick±slip motion, associated to the direction tangential to the motion, is equal to 0 and the second exponent tends to ÿ1. Another characteristic aspect of this system, which highlights its non-smoothness and its fast convergence to the limit cycle, is given by its divergence: Eq. (5) can be written as a system of ®rst order di€erential equations, in the case of …X_ ÿ Vdr † < 0 Eq. (5) becomes: X_ ˆ Y ;

Y_ ˆ ÿX ‡ Md …Y ÿ Vdr †

…7†

with an equivalent relation derivable if …X_ ÿ Vdr † > 0. During the slip phase the divergence is given by div…X_ ; Y_ † ˆ

oX_ oY_ oMd ‡ ˆ oY oX oY

…8†

which is always positive (see Fig. 7) if b ˆ 0 or b small and |Y| not too high. In this way a small area in the phase space is enlarged under the action of the ¯ow. The mathematical model of the friction law concentrates the dissipation on the Y ˆ Vdr line at the end of the slip phase where the quadrilateral area of Fig. 7 is ¯attened onto an interval. It is apparent that these non-smooth systems with variable structure cannot be analysed by methods that require certain smoothness assumptions on the non-linear functions involved [17,18]. For example, the common methods for calculating Lyapunov exponents [19,20] cannot be directly applied: new means of investigation have to be found.

Fig. 7. In the region of the phase space where X_ < Vdr the divergence of the dynamic system is positive, the square area of initial conditions A0 B0 C0 D0 is mapped by the ¯ow to a larger area A1 B1 C1 D1 which is then ¯attened to a segment A2 B2 C2 D2 .

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

297

3. The multi-degree of freedom mechanical system The conventional mechanical model under investigation is composed of a number of blocks arranged along a line (Fig. 8). Three linear springs connect each block to a ®xed body and to the two adjacent blocks while a moving belt supports them. The reference con®guration is chosen as the con®guration for which all the springs assume their natural length and no elastic force is exerted on the blocks. In the reference con®guration . . . Xh ˆ Xi ˆ Xj ˆ


ˆ 0, where Xi indicates the dimensionless displacement of the ith block from the reference con®guration. When the intensities of inertial, friction and elastic forces are such that the blocks ride on the belt with no relative motion between them and the belt (global stick phase), the motion equations for the ith block are: Xi ˆ 0;

X_ i ˆ Vdr ;

Xi ˆ Xi0 ‡ Vdr s;

…9†

Xi0 is the initial value of the displacement of the ith block and s is the dimensionless time variable. This motion augments the elastic forces acting between the ®xed body and the blocks, whereas the length of the coupling springs remains unchanged, up to the point where (generically) one block starts slipping. The geometrical locus, where the slipping condition is veri®ed for the ith block, is given by the equation: ai Xi ‡ ahi …Xi ÿ Xh † ‡ aij …Xi ÿ Xj † ˆ Mst;i ;

…10†

where Mst;i is the ratio between the maximum static friction force acting on the ith block and a reference friction force, ai is the dimensionless stiffness of the spring connecting the ith block to the ®xed support and aij is the dimensionless stiffness of the coupling spring between the ith and the jth block (see Fig. 8). Given a mechanical system with n blocks, at the time that one of Eq. (10) reaches equality, the corresponding block begins to slip with a velocity algebraically lower than the driving velocity. The slipping mode dynamic equation is: Xi ‡ ai Xi ‡ ahi …Xi ÿ Xh † ‡ aij …Xi ÿ Xj † ˆ Md;i …X_ i ÿ Vdr †:

…11†

In Eq. (11) no mass coecient appears because we assume that all blocks have the same mass which, by means of the adimensionalisation procedure, becomes equal to 1 in the dimensionless equations. The slip motion terminates when the velocity of the block is again equal to Vdr ; at this moment two further possibilities are given: (a) if the following condition is true, then a new stick phase begins, ÿMst;i < ai Xi ‡ ahi …Xi ÿ Xh † ‡ aij …Xi ÿ Xj †;

Fig. 8. Multi-block friction oscillator.

…12†

298

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

(b) or, if the following condition is true ÿMst;i > ai Xi ‡ ahi …Xi ÿ Xh † ‡ aij …Xi ÿ Xj †

…13†

then a new slipping phase begins, characterised by the fact that X_ i > Vdr and governed by the equation Xi ‡ ai Xi ‡ ahi …Xi ÿ Xh † ‡ aij …Xi ÿ Xj † ˆ ÿMd;i …jX_ i ÿ Vdr j†:

…14†

The preceding case (a) corresponds to the one shown in Fig. 3 of the one-block system whereas Eq. (13) corresponds to Fig. 4. In any case, if the driving velocity is not too high, the block will end up undergoing a new stick motion. Finally, provided that the driving velocity is low enough, all blocks will stick to the belt at some time and a new global stick phase will once again occur. The graphic representations of Section 5 require the de®nition of a conventional Poincare section. In general the ¯ow of the dynamic stick±slip system is transverse in its 2n-dimensional phase space to the (2n ÿ 1)-dimensional subspace given by X_i ˆ 0. This property immediately suggests that this (2n ÿ 1)-dimensional subspace may be chosen as the Poincare section. Given the derivative: dXi dXi =dt X_ i ˆ ˆ ˆ0 _ _ dXi dXi =dt ÿai Xi ÿ aki …Xi ÿ Xk † ÿ aij …Xi ÿ Xj † ‡ Md;i …X_ i ÿ Vdr †

…15†

it is apparent that the function Xi …X_ i † is perpendicular to the subspace X_ i ˆ 0 satisfying, in general, the transversality conditions required by the de®nition of a Poincare section. The numerical integration of the equations of motion con®rms that a Poincare map can be de®ned by choosing the points where an increasing velocity of the ith block becomes equal to 0. At those points X_i ˆ 0 while the other (2n ÿ 1) lagrangian variables …X1 ; X_ 1 ; . . . ; Xi ; . . . ; Xn X_ n †P constitute the co-ordinates of a (2n ÿ 1)-dimensional Poincare map. In the examples presented in Section 5 the Poincare section is de®ned by the equation X_1 ˆ 0. 4. n ÿ 1 Dimensional maps Intuition suggests, and numerical simulations con®rm, that if the driving velocity is ÔlowÕ, all existing attractors are characterised by global stick phases (from now on g.s.p.) that is by phases of motion during which all the n blocks of the system are sticking on the belt. Therefore we de®ne a driving velocity as low, for a system of n blocks, if all existing attractors are characterised by g.s.p and we limit our study to the system with low driving velocities. The dynamics of this kind of system generates in a simple way a (n ÿ 1)-dimensional discrete map. In fact during a g.s.p. it is clearly possible to compute the following variables: di ˆ Xi‡1 ÿ Xi

…16†

which are constant during the g.s.p. and represent the relative displacement between block i ‡ 1 and block i. After a slip phase, during a new g.s.p., some or all of the variables di assume a value which is different from the previous g.s.p. Therefore the dynamics of the system generates a sequence of successive (n ÿ 1)-ple of values of the following type: . . . …d1 ; d2 ; . . . ; dnÿ1 †k ; …d1 ; d2 ; . . . ; dnÿ1 †k‡1 . . .

…17†

which can be interpreted as the terms of a (n ÿ 1)-dimensional map: d1;k‡1 ˆ f1 …d1 ; d2 ; . . . ; dnÿ1 †k ; d2;k‡1 ˆ f2 …d1 ; d2 ; . . . ; dnÿ1 †k ; ...

…18†

dnÿ1;k‡1 ˆ fnÿ1 …d1 ; d2 ; . . . ; dnÿ1 †k : The (n ÿ 1)-dimensional map is de®ned in a ®nite portion of Rnÿ1 as shown in the following examples. The use of the (n ÿ 1)-dimensional map is clearly more convenient than the use of the common Poincare map

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

299

introduced in the previous section, which is de®ned in a (2n ÿ 1)-dimensional space. Nevertheless the (n ÿ 1)-dimensional map presents some discontinuities, as will be shown in Section 5, so that the dynamics of the mechanical system is more clearly illustrated by the Poincare map. This is the reason why this map is used in some ®gures of Section 5. n ˆ 2. In the case of a two block system [9] it is possible to de®ne a 1-d map: dk‡1 ˆ f(dk ) where d ˆ X2 ÿX1 . The portion of the phase space in which the g.s.p. may occur is given by the relations: ÿ Mst;1 6 a1 X1 ‡ a12 …X1 ÿ X2 †  Mst;1 ; ÿ Mst;2 6 a2 X2 ‡ a12 …X2 ÿ X1 †  Mst;2 :

…19†

This region can be represented in the plane (X1 , X2 ) as a parallelogram (see Fig. 9) and it is the extension to the two-dimensional case of the segment ÿ1 < X < 1 of the single block system. The g.s.p. trajectories are straight lines parallel to the bisection line of the ®rst quadrant. The extreme values of the variable d are assumed at the points A and C of the perimeter of the parallelogram: therefore the map is de®ned in the interval dA < d < dC . Once the de®nition ®eld of the map is known, it is possible to ®nd its form numerically: a dense grid of initial conditions is chosen along the line AC, with initial velocity equal to the driving velocity, and the ®rst iterate of the map is computed for each of them. An example of such a computation is shown in Fig. 10. The map has two periodic attractors, one is a ®xed point, indicated by a circle, the other a period ®ve attractor, indicated by ®ve squares. The intervals along the horizontal axis indicated by the letter a are the basin of attraction of the ®xed point. The remaining parts of the interval dA < d < dC are the basin of attraction of the period ®ve attractor. The map illustrates some important features of the system. Non-smoothness of the forces acting on the blocks and discontinuous de®nition of the motion equations generate a non-continuous map, but the discontinuities are well localised and the map can in fact be assumed piece-wise continuous. Therefore it is reasonable to assume that convergence to the attractors takes place in an in®nite time so that it is possible to compute the Lyapunov exponents. The phase space of a two-block system is four-dimensional, assuming that four Lyapunov exponents can be de®ned; two of them tend to ÿ1 because of the degeneration of motion in a two-dimensional phase space during the g.s.p. A third exponent is zero, in the direction

Fig. 9. Locus of possible stick phases for the two block system.

300

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

Fig. 10. One-dimensional map of the two-block system with parameter values given by: a1 ˆ a2 ˆ 1:0; a12 ˆ 1:2; Mst;1 ˆ 1:0; b ˆ Mst;2 =Mst;2 ˆ 1:3; ) dC ˆ ÿdA ˆ 0:6764 and driving velocity Vdr ˆ 0:12, friction parameters a ˆ 0.0, b ˆ 0.0, c ˆ 3.0.

tangential to the motion, the fourth exponent quali®es the nature of the motion and is assumed to be the exponent of the map. Since the form of the map is not explicitly known but numerically computed, its Lyapunov exponent is also numerically computed by means of a ®nite di€erence technique [1] as: kˆ

m 1X df;i‡1 ÿ dp;i‡1 ; m iˆ1 df;i ÿ dp;i

…20†

where df is the ®ducial trajectory [20], dp a perturbed trajectory close to the ®ducial one and m is a ÔlargeÕ integer number which numerically approximates a limit towards in®nity [21]. In the numerical computations dp;i is chosen in such a way that jdf;i ÿ dp;i j ˆ 10ÿ8 , whereas the numerator of Eq. (20) is computed by integrating the dynamics of the system. n ˆ 3. In the case of three blocks the locus of possible stick phases is given by: ÿ Mst;1  a1 X1 ‡ a12 …X1 ÿ X2 † ‡ a13 …X1 ÿ X3 † 6 Mst;1 ; ÿ Mst;2 6 a2 X2 ‡ a12 …X2 ÿ X1 † ‡ a23 …X2 ÿ X3 † 6 Mst;2 ;

…21†

ÿ Mst;3 6 a3 X3 ‡ a13 …X3 ÿ X1 † ‡ a23 …X3 ÿ X2 † 6 Mst;3 ; the projection of this prism of the three-dimensional space (X1 , X2 , X3 ) on the space of the two variables d 1 , d2 : d1 ˆ X2 ÿ X1 ; d2 ˆ X3 ÿ X1 ;

…22†

gives rise to the hexagon of Fig. 11 which is the ®eld where the two-dimensional map of Eq. (23) is de®ned. d1;k‡1 ˆ f1 …d1 ; d2 †k ; d2;k‡1 ˆ f2 …d1 ; d2 †k :

…23†

Note that the hexagon of this ®gure belongs to the plane of the variables d1 and d2 and therefore corresponds to the segment AC of Fig. 9, whereas the parallelogram of Fig. 9 would correspond to a

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

301

Fig. 11. Field of de®nition of the two-dimensional map generated by the dynamics of a three-block system. The system parameters are: a1 ˆ a2 ˆ a3 ˆ 1:0; a12 ˆ a13 ˆ a23 ˆ 1:0; Mst;1 ˆ 1:0; b1 ˆ Mst;2 =Mst;1 ˆ 1:10; b2 ˆ Mst;3 =Mst;1 ˆ 1:18.

three-dimensional prism in the space of the dof X1 , X2 , X3 . In this case it is not possible to visualise the shape of the map but, assuming that the map possesses the same features as the one-dimensional case, it can be used to ®nd ÔallÕ the attractors, their basins, their two signi®cant Lyapunov exponents and it can also be used to implement some path-following techniques as shown in the next section.

5. Results In this section some results are presented in order to describe the dynamics of multiple friction oscillators. They have been obtained by using the computational techniques described in Ref. [1]. All the attractors of a three block system can be found, in the case of low driving velocity, by investigating the dynamics of the two-dimensional map de®ned in the hexagonal ®eld of Fig. 11. This can be easily done with a cell-mapping method [22] which ®nds all attractors and their basins, as shown in Fig. 12. The cell-mapping technique ®nds two coexisting attractors which are shown in Fig. 13. 800.000 initial conditions are chosen in the ®eld where the 2-D map of Eq. (23) is de®ned. All these motions are integrated in time until they converge onto one of the two coexisting attractors. Initial conditions of motions converging onto the attractor of Fig. 13(a) are indicated by white points in Fig. 12, whereas black points indicate initial conditions of motions converging onto the attractor of Fig. 13(b). Once all the attractors for a given parameter set have been found, the study of the dynamics of the system is usually continued with the investigation of the evolution of such attractors as a parameter is varied. The evolution of the two attractors as the driving velocity changes can be easily followed by means of a continuation routine: the converged values of the phase variables for the previous value of the driving parameter are chosen as initial conditions for a slightly varied value of the same parameter. Figs. 14(a) and 15(a) have been obtained in this way, which is only suitable for following stable motions.The unstable motions can be followed using path-following routines which are often developed for 2-D Poincare maps of three-dimensional non-autonomous smooth systems. These path-following routines can be adapted to follow the two-dimensional map previously de®ned by Eq. (23). In this way it is possible to obtain some

302

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

Fig. 12. Basins of attraction, for the system of Fig. 11, Vdr ˆ 0.22, friction parameters a ˆ 0.0, b ˆ 0.0, c ˆ 3.0.

Fig. 13. Coexisting attractors, parameters as in Fig. 12. (a) Attractor of the white basin of Fig. 12. (b) Attractor of the black basin of Fig. 12.

bifurcation diagrams, as shown in Figs. 14(b) and 15(b) for the di€erent attractors of Fig. 13. In these ®gures the X1P variable is plotted to compare with the continuation diagrams, but actually the two-dimensional map is followed. Finally the spectra of the Lyapunov exponents of the map are easily computed and show the existence of periodic, quasi-periodic, chaotic and hyper-chaotic motions (see Figs. 14(c) and 15(c)). We observe that a motion is de®ned as ÔhyperchaoticÕ if it is characterised by more than one positive Lyapunov exponent [23]. An example of a bifurcation diagram of the map is given in Fig. 16 where the discontinuous nature of the map is apparent. By choosing the periodic motion indicated by the letter a in Fig. 16 (which is motion a of Fig. 13) as the initial motion of the bifurcation path, the corresponding map possesses an

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

303

Fig. 14. Attractor (a) of Fig. 13: (a) continuation diagram, (b) bifurcation diagram, (c) Lyapunov exponents. In ®gure (b), for reasons of clarity the stable motions are not distinguished from the unstable. Such a distinction is made in Fig. 16. Fig. 15. Attractor (b) of Fig. 13: (a) continuation diagram, (b) bifurcation diagram, thick lines indicate attractors, thin lines unstable motions, (c) Lyapunov exponents.

attracting ®xed point. As Vdr is increased the ®xed point becomes unstable at a fold bifurcation for Vdr ˆ 0.223902. The driving velocity of the unstable motion decreases up to a new fold bifurcation at Vdr ˆ 0.178120 which is closely followed by a ¯ip bifurcation at Vdr ˆ 0.178381. The bifurcation path has then a discontinuity at Vdr ˆ 0.181092 (points c) where a new branch of the map appears with the value d2 ˆ 0.306718. For higher values of driving velocity the unstable two-point map exists up to Vdr ˆ 0.210081 (points b) where a new discontinuity interrupts the upper branch of the map. The unstable lower branch exists for higher values of Vdr up to Vdr ˆ 0.218754 where a reverse ¯ip bifurcation takes place and generates a stable one-point map in such a way that the path returns to the starting point a. Figs. 14±16 show that, in spite of being non-smooth, this system seems to be a€ected by standard bifurcations. Moreover the two-dimensional map gives an explanation of the sudden disappearance of the chaotic attractor of Fig. 14(a) for Vdr  0.1625 As shown in Fig. 17 the chaotic attractor collides with basin

304

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

Fig. 16. Attractor (a) of Fig. 13, bifurcation diagram of the 2-D map. 1-point or 2-point map corresponding to the 5-point Poincare map shown in Fig. 14(b), thick lines indicate attractors, thin lines unstable motions. This ®gure shows a Ôclosed loopÕ for the twodimensional map which is not continuous.

Fig. 17. Collision between the chaotic attractor of Fig. 14(a), Vdr  0.1625, and its basin boundary.

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

305

boundaries and this collision is the probable cause of its death since it is believed that a similar collision takes place in the six-dimensional phase space between the attractor and a stable manifold [24]. The attractor is composed of three branches contained in closed curves and indicated by the numbers 1, 2, 3, its basin is the white basin. Comparison with Fig. 12, and many intermediate ®gures which are not shown here for brevity, seem to reveal that no signi®cant basin erosion [25] takes place, when the parameter Vdr changes, before the crisis, as has already been observed in a simpler system [11]. It is worth stressing that all computations are carried out with no simplifying assumption of friction characteristic and motion equations: no smoothing hypothesis [7,26] is introduced but the fully discontinuous stick±slip system is numerically integrated.

6. Conclusions Some features of the dynamics of multiple friction oscillators are described by means of the introduction of ®nite-dimensional discrete maps. In this way it is possible to de®ne a limited region in which all attractors can be located and their basins found. The same maps are used to compute the signi®cant Lyapunov exponents and to implement path following techniques which allow for a wide description of the system dynamics. The dynamics of the system is characterised by periodic, quasi-periodic, chaotic and hyper-chaotic attractors and by standard bifurcation as fold, ¯ip and collision of chaotic attractors with basin boundaries.

References [1] U. Galvanetto, S.R. Bishop, Computational techniques for non-linear dynamics in multiple friction oscillators, Comp. Meth. Appl. Mech. Eng. 163 (1998) 373±382. [2] R. Burridge, L. Knopo€, Model and theoretical seismicity, Bull. Seism. Soc. America 57 (1967) 341±371. [3] S.C. Cohen, Computer simulation of earthquakes, J. Geoph. Res. 82 (1977) 3781±3796. [4] J. Nussbaum, A. Ruina, A two degree-of-freedom earthquake model with static/dynamic friction, PAGEOPH 125 (1987) 629±656. [5] J. Huang, D.L. Turcotte, Evidence for chaotic fault interactions in the seismicity of the San Andreas fault and Nankai trough, Nature 348 (1990) 234±236. [6] K. Popp, P. Stelter, Stick±slip vibrations and chaos, Phil. Trans. Roy. Soc. Lond. A 332 (1990) 89±105. [7] M. Oestreich, N. Hinrichs, K. Popp, Bifurcation and stability analysis for a non-smooth friction oscillator, Arch. Appl. Mech. 66 (1996) 301±314. [8] N. Hinrichs, M. Oestreich, K. Popp, Dynamics of oscillators with impact and friction, Chaos Solitons and Fractals 8 (1997) 535± 558. [9] U. Galvanetto, S.R. Bishop, L. Briseghella, Mechanical stick±slip vibrations, Int. J. Bif. Chaos 5 (1995) 637±651. [10] M. De, Sousa Vieira, Chaos in a simple spring-block system, Phys. Lett. A 198 (1995) 407±414. [11] U. Galvanetto, Bifurcations and chaos in a four-dimensional stick±slip system, J. Sound and Vibration 204 (1997) 690±695. [12] J.J. Stoker, Nonlinear Vibrations, Interscience Publishers, New York, 1950. [13] A.H. Nayfeh, D.T. Mook, Nonlinear Oscillations, Wiley, New York, 1979. [14] J. Argyris, H.-P. Mlejnek, Texts on Computational Mechanics, Vol. V, Dynamics of Structures, North-Holland, Amsterdam, 1991. [15] P. Conti, Sulla resistenza di attrito, Rd. Acc. Lincei 11 (1875). [16] R.H. Abraham, C.D. Shaw, Dynamics ± The geometry of behaviour, Aerial Press, Santa Cruz, 1983±1988. [17] J. Guckenheimer, P. Holmes, Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields, Springer, New York, 1986. [18] J. Argyris, G. Faust, M. Haase, Texts on Computational Mechanics, An Exploration of Chaos, Vol. VII, North-Holland, Amsterdam, 1994. [19] G. Bennettin, L. Galgani, A. Giorgilli, J.M. Strelcyn, Lyapunov characteristic exponents for smooth dynamical systems and for Hamiltonian systems; a method for computing all of them, Meccanica 15 (1980) 9±20. [20] A. Wolf, J.B. Swift, H.L. Swinney, J.A. Vastano, Determining Lyapunov exponents from a time series, Physica D 16 (1985) 285± 317. [21] J.C. Roux, R.H. Simoyiand, H.L. Swinney, Observation of a strange attractor, Physica D 8 (1983) 257±266. [22] T.S. Parkerm, L.O. Chua, Practical Numerical Algorithms for Chaotic Systems, Springer, New York, 1989. [23] O.E. Roessler, An equation for hyperchaos, Phys. Lett. A 71 (1979) 155±158.

306

U. Galvanetto / Comput. Methods Appl. Mech. Engrg. 178 (1999) 291±306

[24] C. Grebogi, E. Ott, Y.A. Yorke, Crises, sudden changes in chaotic attractors and transient chaos, Physica D 7 (1983) 181± 200. [25] M.S. Soliman, J.M.T. Thompson, Integrity measures quantifying the erosion of smooth and fractals basins of attraction, Journal of Sound and Vibration 135 (1989) 453±475. [26] P. Stelter, W. Sextro, Bifurcations in dynamic systems with dry friction, Int. Series Num. Math. 97 (1991) 343±347.