Intermittency in impact oscillators close to resonance

We examine the change in behaviour of the solutions of a simple one degree of freedom, periodically forced, impact oscillator following a grazing bifurcation in which an impact of zero velocity occurs following a change in one of the parameters of the system. It is shown that such a bifurcation leads to intermittent chaotic behaviour with low velocity impacts followed by an irregular sequence of high velocity impacts. We also show that there is a natural, discontinuous one-dimensional map associated with this relating one low velocity impact to the next and the properties of this map are analysed. We also construct the bifurcation diagram of the change in behaviour and show that this contains a series of periodic windows, with the period of the solutions increasing monotonically by one in each successive window as the bifurcation point is approached. By restricting our attention to the resonant case where the forcing frequency is twice the natural frequency of the oscillator it is possible to make asymptotic estimates of the form of the intermittent chaotic behaviour and these estimates are compared with some numerical calculations.


Introduction
In many forced vibrating mechanical systems, oscillating components collide with each other or with rigid obstacles. The motion of the components is then a combination of a smooth motion governed by a differential equation interrupted by a series of non-smooth collisions. Such systems are usually termed impact oscillators and have recently been the subject of growing interest both in the engineering and mathematical literature. The reasons for this are twofold: firstly, they can be used to model a wide variety of different systems arising in engineering, for examples see Thompson and Stewart (1986), Shaw and Holmes (1983), Shaw (198S), and Moon (1987). Secondly, they are examples of nonlinear dynamical systems with discontinuities and as well as exhibiting the usual complexity of behaviour associated with smooth nonlinear systems we also observe new forms of behaviour related to the discontinuities. Some of this is reported in the papers by Whiston (1987), Hindmarsh and Jeffries (1984), Nordmark (1991), Foale and Bishop (1992), Nusse et a1 (1993), Budd et a1 (1993) and Dux (1993).
In this paper we investigate a particularly interesting phenomenon which can occur in an impact oscillator, namely a bifurcation from a steady periodic solution to an intermittent chaotic one as a parameter-in this case the mean separation 0951-7715/94/041191+34$07.50 0 1994 IOP Publishing Ltd and LMS Publishing Ltd between the oscillating components-is vaned. This bifurcation has many similarities to the intermittency phenomena described by Pomeau and Manneville (1980), which arise in smooth systems when a stable fixed point disappean at a saddle-node bifurcation. The presence of the discontinuity in impact oscillators gives rise, however, to new and fascinating behaviour and a new form of bifurcation called a 'grazing' bifurcation. Grazing bifurcations have been identified and studied by Whiston (1987) and Nordmark (1991) for the respective cases of unstable periodic motions corresponding to a saddle point of an associated PoincarB map and stable periodic motions corresponding to a node of the same map. In this paper we extend the discussion of Nordmark (1991) by looking at a special case of a grazing bifurcation in which we can reduce the study of an impact osdlator to that of a one-dimensional map.
A simple model for a forced impact oscillator is a linear one-degree-of-freedom system in which an oscillating component at a position x ( f ) collides with an obstacle at position U. Between impacts the x ( t ) satisfies the simple differential equation (1.2) An impact occurs when x = U and we model this by an instantaneous restitution law: In Nordmark (1991) and Foale and Bishop (1992) somewhat more general systems are studied, but the simplification we make in looking at an undamped piecewise linear system does not significantly detract from the generality of our results. The system above is dissipative as r < l and it is this case that we study in this paper. (The problem with r = 1 is closely related to work on billiaxd balls, see for example the review by Katok and Strelcyn (1986).) Indeed, for all computations and numerical examples we consider the special case r = 0.8. We presume that the system starts in some (bounded) initial state and that we are observing its subsequent motion. It has been shown by Whiston (1987) that when r < 1, there is a bounded set which attracts trajectories for all initial data and we study the motion in this set. A detailed analysis of all subsequent forms of motion and of some of their bifurcations is presented in Budd et al (1993). There, the system is defined to be 'resonant' if o = 2n, (this leads to the largest amplitude motions when U = 0) and we study the case o = 2 in this paper.
The system (Ll), (1.2) leads to a discontinuous dynamical system for two reasons. Firstly, and mwt obviously, there is a discontinuity in the velocity of x at each collision. Secondly, the system has a discontinuous dependence upon the initial data which occurs when it just misses an impact or has an impact with a very low velocity. We can see this by considering figure 1.1 in which three adjacent initial states lead to a motion (A) with a low velocity impact, (B) with a zero velocity impact and (C) with no impact. Although these three trajectories look very similar,  the addition of an impact in (B), (C) profoundly changes the nature of the dynamics and makes the whole system much more sensitive to changes in the initial data. Full details of the mechanism behind this are given in Whiston (1987). This leads to sudden destabilisation of periodic states and the occurrence of chaotic motion. The zero velocity orbit in (B) is termed a graze in Whiston (1987), Nordmark (1991), Foale and Bishop (1992) and the occurrence of a graze as a parameter is varied is termed a grazing bifurcation (in Nusse et a1 1993) a grazing bifurcation is presented as an example of a border-collision).
Because of the dissipation in (l.l), (1.2) the impact oscillator settles down onto an attractor after a transient period of several impacts. The simplest form of such an attractor is a periodic motion in which the system exactly repeats itself between impacts. (These forms of motion are analysed in detail in Budd et a1 (1993).) In figure 1.2 we present an example of a globally attracting, stable periodic orbit which arises when $ = 2, U = -0.33 and r = 0.8. We also show the transient motion of about 30 impacts which arises from taking as initial data, x = U , i = 0 at t = 0, and observe that the periodic motion has a local maximum value between impacts. If U is decreased from -0.33 and w is kept fixed at 2, then we may continuously deform the periodic orbit until this maximum value grazes the obstacle. When r = 0.8, this occurs when U = o* = -0.3312. . . . For U < U* the local maximum is replaced by a low velocity impact and the periodic orbit disappears. For U < U* we observe instead, an intermittent chaotic motion, which for U close to U* remains in a neighbourhood of the original periodic motion. This transition to chaotic behaviour through the grazing of a stable periodic orbit was first observed by Nordmark (1991). The form that this motion takes when o = -0.333 is presented in figure 1.3.
In this motion we see transient behaviour of several high velocity impacts (very similar to that observed for U = -0.33 before the system settled down to a stable periodic state), followed by a low velocity impact, followed by another transient periodic of high velocity impacts and so on. Indeed, we observe N, high velocity impacts between the ith and (i + 1)th low velocity impacts, where the value of Ni depends chaotically upon i, but has a maximum, N ( u ) , of N(u) = 25 impacts. We show in section 4 that and further, that if U((.*) is the impact velocity of the periodic orbit and vmax(u) the maximum impact velocity of the chaotic motion, then where K depends slowly upon U. Result (1.4) was first observed by Nordmark In figure 1.4 we plot one succesive low-velocity impact velocity U,, against the next low velocity impact for the case of U = -0.333. (Here an impact velocity vi is defined to be 'low' if u,<O.5.) This figure is striking as it appears that when U is close to U * , there is a natural, discontinuous, piecewise smooth approximately one-dimensional map effectively relating one low velocity impact to the next. A major purpose of this paper is to construct and analyse this map and show why it is asymptotic to a one-dimensional map in the limit of U+ U*. (In the papers by Foale and Bishop (1992) and Nusse er a1 (1993) an alternative approach to studying grazing is taken in which various one-dimensional maps with discontinuous first derivative are studied and their bifurcations analysed.) In figure 1.5 we present a bifurcation diagram showing the different forms of motion as U is varied, plotting the values of U at each impact. This diagram was computed for each value of U by calculating the first 4000 impacts of the impact oscillator and then recording the velocity of the following 5000 impacts. The most significant feature of this figure is the very rapid increase in the mean velocity of impact as U decreases through U* and the simultaneous appearance low velocity (1991). impacts. Indeed, the velocity increases proportionally to the square-root of ( U * -a ) as predicted by (1.4). This is a significant observation from an engineering point of view, as the rapid increase in the mean velocity will lead to dramatic increase in the mean wear rate of the system. As U becomes more negative we also observe windows of periodic behaviour and period doubling bifurcations. These were also observed and studied for the more general system by Nordmark (1991). This behaviour is reminiscent of that of the one-dimensional logistic map and is partly due to properties of the natural discontinuous one-dimensional map, f, described above.
The layout of this paper is as follows. In section 2, we define a two dimensional (impact) map. P :SI X R+S' X R, associated with (1.1), (1.2) and will determine some of its properties. In particular we study the dynamical effects of grazing and a consequent local reduction of P to a one-dimensional map. (This section is rather technical and some of it could be ommited at a first reading. Further details are given in either Whiston (1987) or Nordmark (1991).) In section 3 we use numerical methods to investigate the bifurcation to intermittent chaos in more detail and demonstrate the existence of an underlying discontinuous one-dimensional map. In this section the main properties of the one-dimensional map are described. In section 4 we examine the local structure of the bifurcation close to U* and prove the results of scction 3, in particular justifying the existence and the properties of the one-dimensional map. Finally in section 5 we extend our analysis to more negative values of U which are not so close to U*. At these values certain global properties of P become important, causing the form of the map to change. This enables us to explain features of figure 1.5 not predicted by the local analysis of section 4.

The impact map P
The basic dynamics of the system (1.1), (1.2) may be studied by reducing it to a map P , which we call the impact map, on the surface {x = a , U > 0) taking the phase and velocity ( + , , . U , , ) at an impact at time t = t,) to the phase and velocity ( + , , U , ) at the next. Here the phase, + is defined by + = tmod(2x,,,,p The geometry of this map has been studied in detail by Whiston (1987) and Chillingworth (1989). Given an impact with phase +,) and velocity U,,, we firstly apply the instantaneous linear impact rule, which can be expressed as a map such that and then allow the system to evolve according to equation (1.1) until the next impact. If F is the map for this evolution then

C Budd and F Dux
There are certain technical difficulties associated with the definition of P when U = 0.
In particular, if U = 0 and the acceleration toward the obstacle is positive then the particle remains stuck to the obstacle for a non-zero period of time. Moreover, if U = O then there are certain values of the phase for which a subsequent impact does not occur and P is not defined. Extensions to the delinition of P which overcome these difficulties are given in Whiston (1987) and their dynamical consequences discussed in Budd and Dux (1994).
The most significant feature of P of interest to us here is the fact, discussed in the introduction, that it is discontinuous. If we define to be the set of zero-velocity or grazing impacts which lie in the range of P and correspond to maxima of x(t), then it is the existence impacts which leads to discontinuities in P.
Lemma 2.1. Let q be a point such that q non-degenerate local maximum at this point, then P is discontinuous at q. I,,, P(q) E I,, and the flow has a Proof. Let q = (CO, U ( , ) where U" # 0 and P(q) = (+,, 0). Furthermore let the flow starting at q be x(t) such that x(+,)) = cr. Then, x(t) has a local maximum at tI (with corresponding phase 4,) such that x ( f , ) = U and Z(fl) = -2m <O. We now consider a sequence of trajectories x(r; (Y) such that dx(t,; a)/dt = 0 and x ( t , ; a) = U + a. If a is close to 0 then, by the continuity of the flow (reversing the direction of time), and the assumption that ul,#O, it follows that x ( f , a ) = u at a time f ( a ) which is a continuous function of a with t(0) = +, > Moreover, the velocity U(.) at the time t(.) will also be a continuous function of a. Such a trajectory will thus correspond to an initial point g(a)=(r(a),u(a)) such that t ( a ) + q as a+O. Now by the continuity of the flow it follows that there is a neighbourhood [Iz, r3] of r , (which does not depend on a if a is small) such that if we take (Y sufficiently small then -3m <ji(r; a) < -m if r2 < r < t3. Thus, if a < 0 then x(r; a) # U at any point in this interval wheras if a>O then x ( r ; a ) = o at a time T ( a ) in the interval t l -(2alm)'" < T < tl -(2a/3m)'". We conclude that if a + 0' then a trajectory starting at the point I(&) impacts at the time Tclose to I , and hence P(f)+P(q) as <(ol)+q. In contrast, if a-tO-, then the trajectory impacts at a time T >r3 which is bounded away from f, as E(a)+q. The map P is therefore discontinuous at the point q. 0 Nofe. We can easily extend this result by defining two disjoint non-empty open subsets VI, V, of a neighbourhood U of q, such that if 5 E V, has a corresponding flow x ( t , f ) with a maximum at a time close to tI with height hm&) then The map P then has a discontinuity as 5 passes from VI to V,.
Corollary 2.2. Let 5 E V, tend to q. If P 2 ( q ) e &, then p(f)-+pZ(q) as <+a Proof. If the height of the maximum close to f, is smaller than but close to U then the proof of lemma 2.1 indicates that the next impact will be at a time T > r, . As the flow is continuous, this impact must be in turn close to a point of impact of the trajectory x ( t ) starting at q, which as P 2 ( q ) must also be the next impact. The result then follows. 0 At first sight it might appear that these discontinuities are unimportant and are merely the bi-product of a 'poor' definition of P. However, a grazing impact profoundly alters the stability of the flow x ( t ) leading directly to the curious dynamics discussed in this paper. Following Whiston (1987). We define the discontinuity set of P to be the set S curve but may contain branch points, which are singularities of the projection. In the present context is the effect of the discontinuity on the dynamics of P, governed by the behaviour of the derivative DP near S which is of interest. More details on the form of S are given in Whiston (1987) and Chillingworth (1989).
Initially we prove a few simple geometric results on S and W.
Lemma 2.3. Where it exists, the map P-l is given by, Reversing the time direction in equation (1.1) or (1.2) has the effect of reflecting the (x, 4, U) phase space about the plane q5 = 0. Thus the map C Budd and F Dux takes the phase and velocity immediately before one impact to the phase and velocity immediately after the previous impact. Pre-multiplying by G-' therefore gives P-'. 0 We deduce that

.
Since S and Ware parameterized by 4 E (@", 4'), we may introduce orientations on them by identifying as the positive direction, the direction in which q5 increases.
With respect to this orientation we consistently define (at least locally) a right-hand side and a left-hand side. As (+o,~o) approaches S from one side then ul tends smoothly to 0 and the image point P ( 4 , , u l ) lies on one side of W. Following Nordmark (1991) we call this side of S (and W) the impact side. Conversly, (+", UO) approaches S from the other side, then + L and uI jump to the other side of W as a new intervening impact occurs. This side of S (and W) is defined to be the non-impacf side. The map P introduces stretching to the phase space on the impact side of S.
For the remainder of this section we consider the behaviour of the map P in a neighbourhood of the set S and in particular consider the stretching introduced by the discontinuity. The discussion in this section will be limited and further details are given in Whiston (1987). We presume that an impact at phase & and velocity vo leads to a second impact at phase 4' and velocity u l . On the impact side of S,ul is close to zero and it is bounded away from zero on the nonimpact side. By solving the differential equations in (1.1) between impacts we may construct an explicit expression €or the Jacobian of P. This i s derived in Budd et d (1993) and has the following form sin(Ao)) -r(cos(&,) --sm(Ao)) AI .

U1
Here A, = cos(04,) -o is the acceleration at each impact and A. = t1to.
Thus, as r < 1, the map P on average is area contracting on phase space. However, if u1 is small, as it will be on the impact side of S , the locally P will introduce considerable stretching.

U
This result may be obtained by a direct calculation but also applies to more general impacting systems where there is no energy loss between impacts. A more general derivation is given in Budd et a1 (1993) where it is also shown that if P is considered to be a map from (+,U*) to (+,U') then Det(Do)=r2. More insight may be obtained by considering a geometrically motivated decomposition of Do motivated by the following result.

It now follows by inspection that
(2.4) Thus as we approach S from the impact side, Do comprises two components of which one becomes singular as ul+O. However, if we approach S from the non-impact side then U , does not tend to zero and Do does not become singular. It is clear from this decomposition of its Jacobian that as ul-tO the map P introduces considerable stretching. Indeed if we consider a point close to S lying on S,, then any vector through this point in a direction not tangent to S, is to leading order stretched by a factor proportional to l / c = l / u , and mapped onto the vector [l,A,]'. Similarly a vector tangent to S, is not stretched and mapped parallel to [l, 01' .
Using this decomposition we may derive an approximate formula for u I as (0)) is a well defined vector normal to S. We now consider a point s on S and a curve of points (+o, u0) on the impact side of S given by (+o, v0) = s -w o ) , ow IMO), p(o))i2.
From the derivation of Do and the formula for (a(c), p(c)) it follows that, Hence, as e + 0 C Budd and F Dux vI(e) = ( 2 e y -+ o(e).

(2.5)
We now consider the effect of applying the map P twice to points on the impact side of S. from the definitions it follows that S is mapped to W under the action of Pz, indeed, points in a neighbourhood of S are strongIy contracted onto W. We see this as follows where MI is a matrix which is continuous as v1 -+ 0.
The following result now follows immediately.
Lemma 2 . 7 . Ifs is a point on Z, then any vector through s is mapped by P onto a tanjent vector of W. Moreover,. the vector b = [-r,AJ' is mapped to zero.
We conclude from this result that a neighbourhood of s is a mapped to a set which is tangent to W. This result will become very important in our discussions of intermittency because the contractive effect of P in the direction of the vector b effectively reduces it to a one-dimensional map in a neighbourhood of s.
We now consider the map Pz. Combining our two expressions for the Jacobian of P it follows that on the impact side of S where the constants a, 6 , c, d are continuous as q-0. side of S by ignoring the collision at D ( P ) would be regular and equal We can use this expression to determine the form of D(P) on the nonimpacting u l ) and setting r = -1 in (2.6). In this case We now use formulae (2.6), (2.7) to determine the image of a rectangular neighbourhood of a point s on S. In particular, we consider a small set K of 'sue' E << 1, parameterized by the scalars 0 and $ such that where (c,,P)' is, as before, the normal vector to S at s, pointing into the impact side. If E > 0 then K lies on the impact side of S and using the previously determined expression for v1 and the above expression for D(Pz) it follows immediately that there is a constant fsuch that where (7, is, as before, the tangent vector to W at P"(s). The rectangle K is thus stretched by a factor and mapped to a long thin set tangent to Wand lying on the impact side of W. Conversely, if E < 0 then K lies on the non-impacting side of S and is mapped by the now regular map to a set of similar dimensions to K on the non-impacting side of W. We illustrate the set K and its images in figure 2.1. We summarize the results of this section in the following theorem.
Theorem 2.8. Let K be the set defined above then K is mapped to a set transverse to W.
vectors normal and tangent to S, then the image of ($,U) under the action of Pz may be expressed in terms of coordinates y and x with respect to the normal and tangent to W so that (1991)  Y I =AB1, yz=A&,x = BE'nB:" + (241 BE'"@^ + C&

AC AC B
The pre-image of this set is to leading order the curve U, given by I"

In
It can be shown easily that B , C>O. The collection of the sets U, as x varies over P 2 ( K ) defines a foliation of K parameterized by x. We define Q as follows Q((0, I))) = (€e*, 0):both (€e*, 0) and e(& $) lie in the same set U, .
Alternatively, to leading order

B
It is clear from its construction that lp2((@, @)I -J'YQ((@, @))) We see from this last result that if E is small then the action of Pz on the impact side of S is essentially equivalent to that of a one-dimensional map. We shall return to this again in section 4 when we study the reduction of P to a one-dimensional map close to the bifurcation.

A qualitative account of intermittency after a grazing bifurcation close to resonance
Having discussed the geometry of the map P we now consider how this geometry determines the dynamics of the impact oscillator. The simplest form of attracting set of the map P is a periodic orbit. Such an orbit is a finite sequence of points, (+(I), U(')), . . . , (+@'I, U("')), such that P(+"', dm)) = (+(I), U(')). To represent a physical solution of (Ll), (1.2) the corresponding trajectory must also satisfy the following two conditions: As the parameters w and c a r e varied, a solution of (3.1) may cease to represent a physical solution of (l.l), (1.2) by violating either one of these conditions in a grazing bifurcation. This occurs in two ways. Either one of the impact velocities becomes negative, violating condition (Ci), or condition (Cii) is violated and an additional impact is included in the sequence.. This latter phenomenon occurs when one element of the sequence (+('I, U(')) crosses the set S. We now study an example of this second form of grazing bifurcation.
It is shown by Budd et a1 1993 that any periodic orbit of the harmonically forced impact oscillator must also be a harmonic or subharmonic of the forcing. Since it must also be a periodic point of P, it will be characterized by a pair of numbers ( m , n), where m is as above and n is the order of the subharmonic of the forcing frequency. In the case of (1,n) orbits, which correspond to fixed points of P, the time between impacts equals n times the forcing period. From this condition we deduce the following formula giving the phase & and the velocity U, of such orbits as (implicit) functions of U and u (see Budd et ai (1993), Whiston (1987), and Hindmarsh and Jeffries (1984) where c, =cos(?).  We parameterize the points on the ellipse in terms of +n, which ranges from 0 at orbit must either intersect the obstacle between impacts or lie completely within the obstacle. Thus, if for some value of U the solution of (3.1) gives a physical orbit, then, by continuity in U it follows that there must be a first point U* > -(yI, with corresponding +,*,u:>O for which condition (Cii) is violated and a grazing bifurcation must occur. For many values of o, this occurs at a value of U* for which the orbit has already lost stability through a period-doubling or a saddle-node bifurcation. However, this is not the case when o is close to the subharmonic value 2n. It was shown by Budd et a1 (1993) that, for o <2n, the major axis of the ellipse has positive slope, for o > 2n, it has negative slope and for w = 2n, it is vertical.
Moreover, when w = 2 n , the value of U, at both the saddle-node and period doubling bifurcations is zero. We refer to the case of w = 2n as the resonanf case.
Indeed, if u=O this leads to the largest amplitude behaviour. These values are significant as we may conclude from the above that: if w is less than, but close to 2n, the fixed point (+,,,u,J will remain stable for u>u*, encounters a grazing bifurcation at U = U* and for U < U* is replaced by a (2,2n) orbit which theorem 2.1 indicates will be unstable.
The numerical studies presented in the introduction showed that for w = 2 and r=0.8, there is a (1,l) orbit which encounters such a grazing bifurcation at U* = -0.331 269 43. . . leading to a chaotic motion when U < U* with a characteristic structure akin to the 'intermittency' phenomenon of Pomeau and Manneville To study the nature of this motion in more detail we record the velocities U, at each impact. In figure 3.2 we present a plot of 1000 values of U; recorded against i (1980).  3 E typical point in the set 2 returns to it under iteration after a number of high velocity impacts, during which it moves toward the position that the k e d point of P occupied in phase space before the grazing bifurcation. The strange attractor has N one-dimensional 'fingers' which intersect close to the earlier position of the fixed point. The effect of P is to map one finger into the next with an individual point 'moving down' the fingers toward the point of intersection from which it is mapped into 2 and then back onto the first finger. Similar figures to this are presented in Nordmark (1991).
In a sense which we shall make precise in section 4, the points in the set 9 are described by a natural coordinate s, so that s = 0 when U = 0 and s increases as U jncreases. In terms of this coordinate we may describe the above map from 2 to itself in terms of the mapf(s). The map F(u) is an approximation to this map due to the projection of dp onto the u-coordinate. By studying the form of the map F(u) indicated by figure 1.4 we may infer various properties of the map f(s) which we list as follows.
1. The map f(s) is a piecewise continuous function. The domain 2 of fis divided into N + 1 subintervals Jk, 0 < k < N , on which f i s continuous, so that Jk = (ak, ak+iI. Here 0 = a , < a, <. . . < aN+l. The index, k , of each interval Jk is such that k + 1 equals the number of iterates of P (or high velocity impacts), which separate the corresponding elements of the set 2. d of length N + 1 or N impacts much more often than shorter excursions. This is illustrated in figure 3.6 which shows a histogram of the numbers of excursions of different lengths for U = -0.333 for which N = 24. 5. In each subinterval Jk, f(s) has a single fixed point corresponding to a ( k + 1, k) periodic orbit of P. If E is small then the slope of f(s) is everywhere less than -1 and thus all the corresponding k e d points are unstable.

6.
As Iu-u*~ is decreased, the number, N + 1 , of subintervals Jk, increases without hound. New intervals appear at a sequence of discrete values U,,, of U.
At each U, a new fixed point appears, corresponding to a grazing bifurcation of a (N + 1, N) orbit of P violating condition (Cii). As v+ v* A schematic illustration of the map f(s) is given in figure 3.7. It should be emphasized that f ( s ) is in fact only an approximation to a one-dimensional map. Indeed, if E =(a*a) then, in the limit of small E, each of the curves defining f(s) will be a set of width O(E~") and length O(e'"). As the ratio of the width to the length is O ( E ) then this set is asymptotic to a one-dimensional set. The onedimensional. nature of the map is also emphasised by the concentration of the measure of the attractor on those impacts corresponding to points in .IN and JN--l.
We note from 5 that on each interval Jk, the mapf(s) has a fixed point which will (in general) be unstable. (The cases where they are stable are discussed in section 5.) Each point corresponds to an unstable (k + 1, k ) periodic solution of the impact oscillator system. It follows from 6 that the grazing bifurcation creates approximately log(~)/(2 log(r)) such points. All of the above results will be proved in section 4 for the case of a close to a*. For more negative values of a, the form of f(s) becomes more complicated than that suggested above and, indeed, loses its 'one-dimensional' character. In figure 3.8 we demonstrate this by plotting the related map F(u) for CT = -0.3375. Although F(u) has broadly similar features to the previous map, the points do not lie on as sharply defined curves but on rather broader sets. Moreover the slope of F(u) on the interval .IN becomes flatter, allowing the fixed point in that interval to become stable for some values of a which then loses stability in a super-critical period doubling bifurcation. This leads to the bifurcation diagram presented in figure 1.5. The reasons for this behaviour are rather subtle and are discussed in section 5.

The form of the grazing bifurcation when U is dose to u*
We now study the system described in section 3 for a close to 'a*, (a*= -0.331 269 43.. . when o = 2 and I = 0.8) to determine a description of the bifurcation to chaos in terms of E = a*a, when E is small. Our analysis is similar to that of Nordmark (1991) but, by restricting ourselves to the near resonant case, we may extend Nordmark's calculations to a detailed investigation of the map f(s).
Our calculations are based upon the observation that if E is small then we may partially analyse the dynamics of the map P in terms of the local linearization of the dynamics of a related map Q in the neighbourhood of a fixed point of this latter map.
We shall assume that a stable (1,l) orbit corresponding to a fixed point y = (+n, U,) of P, undergoes a grazing bifurcation when U = U*. At this value, y crosses a point of intersection P = (+, , U,) E w n s at non-zero speed such that y lies above p (i.e. U,, >U,) if U > U*.
Our discussion involves certain technical conjectures about the generic behaviour of W and S as functions of E , which have been verified numerically for the range of parameters considered. Figure 4.1 shows the curves W and S for w = 2, and for U in a neighbourhood of U*. These curves intersect at a single point, p , and the relative orientations are such that the impact side of S lies below S and to the left wheras the impact side of W lies below Wand to the right. The se& Wand Shave a particularly simple form at these parameter values, with no branch points, and this simplicity makes analysis significantly easier. For other parameter ranges where chaos occurs, Wand S may have a much more complicated geometry and the chaotic attractors are correspondingly harder to analyse; some examples of these are given in Budd et al(1993) and in Whiston (1987).
A central idea in this section is the construction of a map Q which is related to P and which perturbs smoothly as E passes through zero.  (ii) if U E S n U, then

Q(u) = P'(u);
(ii) if U E U and U is on the impact side of S, then where 4' and U' are the phase and velocity of the second non-negative velocity crossing of x = U by the non-impacting solution ~( t ) of (1.1) (in other words, we allow one penetrating excursion before the next impact), Proof. The construction of Q is identical to that of P , except that an impact is ignored. Hence the continuity and differentiability properties are similar provided that Q(u) does not lie in &.
The map Q is illustrated in figure 4.2.
The following results follow immediately from the definition of Q. (iii) For U S c*, the non-physical fixed point Y ( E ) , obtained from solving the algebraic conditions (3.1), (3.2), is a fixed point of Q so long as it lies in U. Moreover, y (~) varies smoothly as E passes through zero and its stability is unchanged.
(iv) If a path r: [0,1] + U crosses S at a point g E U, with g = r(t8), r(0) on the non-impact side and r(1) on the impact side. Then we have that 0 lim P(r(t)) = P'(g) = Q(g). ' Ti , The significance of the map Q is that points on the non-impact side of S behave identically under the iteration of both P and Q. In particular, when U < U * , the fixed point, y , of Q, has a basin of attraction of which part lies on the non-impacting side of S and hence points on this side will be attracted toward y under the iteration of P until they cross S. During this motion, when they are close to y their motion can be accurately described by the local linearization of Q about y . By choosing E small we can always assume that this will be the case.
We now study some of the aspects of the local dynamics of Q in the For w = 2n, it follows from the expression of the Jacobian given in lemma 2.6 neighbourhood of y. that as sin(A,) = 0 and cos(A,) =~-l then the Jacobian of Q at y = (+n, V , ) is This has a single eigenvector ( : ) with a double eigenvalue. At resonance, the pointy is thus a degenerate node with nearby points being attracted along a one-dimensional local stable manifold. For the case considered, u<O and A >O. The dynamics of Q close to y is thus locally approximated by the pair of linear difference equations xk+I = rxh yk+l = (1 + r)AXh + rYk where xk = @h -@n and yk = uh -U,. Solving these gives k xk = r xu Moreover, Q has a one-dimensional invariant manifold, which coincides exactly . We see this as follows.
with the vertical one dimensional eigenspace spanned by

U@,) = I v o -2w-y S i n ( W 4 " )
Thus, for w = 2n, an impact at (&, uu) will be mapped by Q to an impact with the 0 We now presume that for sufficiently small E = U* -U , the set S intersects C at a point (x, y ) = (4 -4n, U -U,,) = (0, ae) where a > 0. Moreover, we shall assume that E is small enough so that this point lies in the neighbourhood of Y(E) for which the local linerization accurately captures the dynamics of Q. We presume further (on the numerical evidence presented in figure 4.1) that the local form of S is given by

(4.3)
We observe that these sets all intersect at the same point when E = 0. We also note, from section 2, that this description of the sets is only valid for those points ( x , y ) which lie above W (i.e. those points which are images under P of points on the nonimpacting side of S . ) We illustrate the resulting sets in figure 4.3.
We now take a similar approach to Nordmark, and construct a trapping region H for the map P by considering sets invariant under P and its iterates. Proof. We presume that (+.U) lies sufficiently close to y so that it is both in its domain of attraction and in a neighbourhood where the local linearisation of Q applies. Thus and hence there must be at least N < m, such that Q"'(4, U) lies on the impact side of S. It further follows from section 2 that, since QN-'(q4,u) lies on the non-impact side of S, QN(4, U ) must lie on the non-impact side (i.e. to the left) of W. Finally, we note that, as (4, U ) lies to the right of C, then the local linearisation of Q implies that e' (+, U ) lies to the right of C for all k > 0. Hence, QN(+, U ) lies to the left of W , the right of C and on the impacting side of S. Consequently it lies in the set E. 0 C Budd and F Dux-

Q k ( 4 , u ) + Y
ask-w, The set E can now be extended to form the trapping region H. 8 -E is small because of the tangency of Pz(C) to Wand hence the point e is given approximately by the intersection between Wand S. It therefore follows from our earlier calculation that (4.4) In particular, if we parameterize W by the scalar t so that W = { ( x , y ) : x = t r b , y = t ( ( l + r ) A b -r c ) , t ER} (4.5) then c i s approximately given when

t =
We note further that b is given when t = 0. Now, the map Pz is a diffeomorphism on 9 and hence maps the boundary of 8 to the boundary of P ' ( 9 ) and the interior of 8 to the interior of P'(9). As 8 intersects S, lies on the impact side of S and is of 'size' E , we may apply theorem 2.8, to conclude that P2 maps 9 to a set of length O(eln) and width O(e3n) approximately tangent to W at b, lying to the right of W and pointing into the non-impact side of S. The set P z ( 9 ) is also illustrated in figure 4.4 and its asymptotic approximation to a one-dimensional set in the limit of € 4 0 leads directly to the form of the map f(s) described in section 3.
Using these results we construct a trapping region for P. Proof. Firstly, we prove that every point of P z ( 9 ) lying on the impact side of S in fact lies in 9. As 8 is bounded by C, S and P'(C), it follows that P'(9) is bounded by Pz(C), Wand PZ(bc) (which is a subset of P"(C)). Hence,Pz(8), must lie in the region bounded by P z ( C ) and W. The intersection of this region with the impact side of S then necessarily lies in 8. We can be more precise by invoking the local linearization for small e Using the approximation (4.4), (4.6) for e and using the local linearisation of Q implies that the point P'(c) approximately lies on the intersection of and Wand is given when Thus, for small e,P2(c) lies strictly between b and c.
Next we consider those points x E P'(9) n S lying on the intersection of P'(9) and S. We know that P(x) E &= {(+, u ) : u = 0}, and that P2(x) = Q(x) E W.

Moreover, x lies between the points p and c on S and thus P2(x) lies between P'(c)
and P'(p) on W. Consequently, we have Pz(x) E 9.
Finally, we consider those points in P2(9) lying on the non-impact side of S. We now presume that E is sufficiently small such that points of O ( d n ) distant from the pointy(€) still lie sufficiently close to y to lie within its domain of attraction under the action of the map Q. It follows from lemma 4.4 that each such point, x, will lie in E c 9 after N ( e ; x ) further iterations. By simple continuity and compactness arguments it follows further that we may take a maximum value of N to be N(E) as x U Thus H is a compact invariant set for P. Lemma 2.4 implies that P will be contracting on Hand, hence, that H contains a set of zero Lebesgue measure, which is the o-limit set of all the points in H.
In  1 S k G N -1.
The sets I, are illustrated in figure 4.6. As E is small they are line segments to a 6rst approximation as Considering the image each set Ik for k > 0, we define Kk c 9 = Ph(lk). As the points in the set 1, lie close to Wit follows that the points in Kk lie close to W"l.
Moreover, for 1 S k S N -1 it follows that Kk is close to 9 f l Wkc' and is bounded above by, and includes, the set Sand is bounded below by, and does not include, the set W. The set KN is close to a subset of 9 n WNC1 and is bounded below by the set W alone. The points where W intersects K, are computable from the local linearization of Q, so that, The sets Kk are also illustrated in figure 4.6.
We now consider the action of P approximately triangular set 9 given by P ( ( -r , A , ) and Al is the acceleration at w. It follows from theorem 2.8 that the vector a is approximately parallel to 9, and that the map P has an approximately zero stretching rate in the direction b and an approximately constant stretching rate in the direction a. Thus we may consider the action of P on 9 to be a map on the single coordinate s. As a first consequence we may estimate the length (and the corresponding value of s) of each set Jk along W2 by rescaling the values o f t corresponding to the intervals Ih. We conclude that there is a constant J such that (4.9) Points in each set . Ik are mapped by P through Ik and D to return to P ( 9 ) in k + 2 iterations. We now consider the behaviour of the iterates of P on each Jk. Each set Jk is bounded below, and does not include, the set Lh = S"' f l P(9). It follows immediately that all points in Lk have the same coordinate s = a h . Similarly, Jk is bounded above, and includes, the set L k + l -S k + 2 n P ( 9 ) for which s = a h t l .
Following k + 1 iterations of P, J, is mapped into the .approximate line segment K,, the points on La lie on K, n W and those points in a neighbourhood of LkCl are mapped to points in a neighbourhood of K k n S . Because the map Pkcl has an approximately uniform stretching rate throughout the set Jh it follows that there is an approximately linear map between Jk and K,. As K, lies on the impact side of S and is always close to S we may now apply the result of theorem 2.8 to deduce. the

Q(sk).)
The set Kk can thus be parameterized in terms of the scalar A so that K,=sk+A (wx-s,) where O S A < 1 if k < N , or A* < A < 1 if k = N for some A* >O.
parameterised by A lies on a line segment interior to 2 such that Invoking theorem 2.8 it now follows that the image of under P of a point P ( s~ + A ( w~ -~k ) ) P(sk) + A"?k (4.10) where P(sk) E Z0 lies on the 'base' of P ( 9 ) and j , = P(wk) -P(sk) is a point near the 'peak' of P ( 9 ) . The sets Jk and the corresponding image sets P(&) are illustrated in figure 4.7.
We may now construct the (approximately) one-dimensional map f ( s ) : 2 -t 2. Each point in 2 has coordinates s and U but because of the zero stretching of P in the direction b the image of such a point is to a very good approximation determined only by the value of s. A typical point in a set Jk has a coordinate s lying in the interval [ak, u k i l ) and is mapped, by Pk into the set Kk and then by P into the set 2 to give a new value of s. This defines the map f(s). Now the map Pk is approximately linear in Jk and hence the value of A corresponding to s is given approximately by

A = ( s -a k ) / ( a , + i -a k ) .
If we combine this approximation with the square-root map in (4.10) we see that f(s) is given approximately by

f ( s ) = . k ( a k + l -s ) ' " / ( a k + l -a k ) ' " , s ( a k ? a k + I ]
where (apart from uo) the values of a, increase with k and are bounded above.
Here we note that on JN the description of f(s) only applies to a subset of the interval (aN, aN+,]. This description of f(s) combined with the explicit calculation of ah from the length of Jh is that given in section 3 and illustrated in figure 3.7.
The approximation of f(s) as a one-dimensional map relies on the result, immediately following from theorem 2.8, that two points p , , p 2 with the same coordinate s and with coordinates U , and u2 with /U, -u2/ = O(e) will be mapped to two points q, and q2 with 19, -q21 = O ( E~~) .
Thus the approximation of f(s) as a one-dimensional map becomes increasingly accurate as E+O. (We note further that it can be shown that [u, -u21 = O((1r )~) and hence this approximation is even more accurate if r is close to 1.) In figure 1.4 we displayed a map F(u):2-tdZ; mapping one low velocity impact to the next. There are difficulties in extending f(s) to such a map because of the non zero width of O(e) of dZ; implies that two points with the same velocity may correspond to different values of s. However, the numerical values given in figure 1.4 strongly imply that this is possible. The numerical observations follow, in fact, from the nonuniform measure on the attractor of the intermittent motion. It is not unreasonable to assume, and the numerical evidence strongly implies, that the distribution of the points within each individual Kk is uniform, that is given that x E &, the probability that the parameter A corresponding to xis less than A is just A. This in turn implies, after applying the map P, that the probability that a point p ( x ) E 2 has velocity less than V is proportional to V2. This again agrees with the numerical evidence. Hence, by far the greatest number of points will be concentrated in the sets Jh where k is close to N and these will be mapped to the line segments P(K,) where k is also close to N. It follows from our earlier description of Kk that for N large these line segments are very close to P(C). The figure 1.4 obtained by plotting points on the attractor is therefore essentially just the map of the single line semgent P(C) into 9 and is necessarily one-dimensional. (We can also see from this figure that the points are concentrated into the sets Jk for large k.)

Global features of the bifurcation diagram
In this section, we study in more detail the features of the bifurcation diagrams presented in figure 1.5 which can effectively be regarded as a bifurcation diagram for the map f(s) when r = 0.8. The discussions in this section will be necessarily more approximate and heuristic than those presented in section 4 as we consider the dynamics of the impact oscillator for larger values of E.
There are three features of figure 1.5 which cannot be explained by the local analysis developed in section 4, and which depend on global features of the impact oscillator which have not yet been considered. These features as follows.
(i) The maximum impact velocity v::,~ of the high velocity impacts in the chaotic motion grows like E" for small values of E , in accordance with the local theory, but at e = 0.003 (for r = 0.8), u:.,~ is asymptotic to an approximate value of The maximum impact velocity, U " , . . , of the low velocity impacts, continues to grow (ii) There are a sequence of windows in the bifurcation diagram in which periodic rather than chaotic behaviour is observed. As E tends to zero, the number of impacts observed in the sucessive windows increases monotonically by one reaching a maximum finite number N,,,. For E close to zero, only chaos is observed. If e* is the greatest value of E such that a k-impact orbit is observed then E k + l / E k % r.
(iii) In each window there is a short period doubling cascade in which the k impact orbit loses stability as E is decreased to a sucession of 2'k-impact orbits. However, after a finite number of these bifurcations the stable 2'k-impact orbit is destroyed in a saddle-node bifurcation and the system becomes chaotic. An example of this is given in figure 5.2 which shows a close up of the bifurcation diagram for (iv) In figure 5.3 we present a bifurcation diagram for a grazing bifurcation when r = 0.3 in which we see a similar picture to that in figure 1.5 but in which the windows of periodic motion are larger and also, although there is a chaotic region for quite small values of E , the motion ultimately stays periodic as ~-0 .
We can explain these features heuristically in terms of the features of the one-dimensional map f ( s ) and global features of the flow. (An alternative explanation is given by Nusse et a1 1993 in terms of a single map with discontinuous first derivative.) As we have seen, when E is small, the map f(s) comprises a sequence of approximately parabolic curves on the intervals Jk = (ak, ak+J, each of which intersects the line f ( s ) = s. These fixed points correspond to a (k + 1, k ) periodic orbit of the impact oscillator. The modulus of the slope of f ( s ) at the points of intersection decreases as k increases and is bounded below by lf'(aN)I giving an -0.3372 < U < -0.3367. will always be unstable if r > $. This explains why we observe an immediate bifurcation into chaos if r = 0.8 but not so if r = 0.3. This accounts for observation (iv).

C Budd and F Dux
The computations illustrated in figures 1.5 and 5.1 show, in contrast to the above, that for n S N,,, = 21, there are windows of values E,* of E , E,* E (en, such that for E E ( E -, E,*) the fixed points are stable. The explanation for this behaviour lies in the change of the global behaviour of the set Pz(D) as E is increased. Providing that E is not too large, the set D continues to be approximately triangular with sides of length proportional to E. Moreover, the set 2 continues to approximate a subset of W 2 and has length proportional to This is in accord with the numerical observations recorded in figure 5.1 that utax is proportinal to E' ' for a wide range of values of E . The set PQ) then continues to closely approximate a subset W3. For very small E this subset of W' is approximately a straight line segment nearly tangent to W. In figure 5.4 we show the global form of the set W3 in which is initially tangent to W but which eventually curves round to intersect W again transversely. The maximum value of U taken by points lying in W 3 is approximately 1.46. This maximum point is in turn mapped by P to a local maximum point on each of the sets W*. When k = 4 or 5 this point is in fact a global maximum value for U during the iteration of P and this global maximum is equal to uk,,. This effectively sets a limit on the maximum value of U on P z ( D ) and hence also on H , providing an explanation for observation (i).
To explain the appearance of the periodic orbits discussed in observation (ii), we consider the effect of the global geometry of the set W 3 on the function f(s), assuming that E is now sufficiently large for Pz(D) to be mapped onto the curved section of this set. We presume that the set W 3 is intersected by the sets SI,. . . , SN, where N decreases as E increases, and that the curvature of W 3 implies that it wiU intersect the set SN (and possibly the sets Sj €or j close to N ) more than once. If E is not too small, then the set J will be sufficiently large for part of it to be mapped under P onto the curved part of W3. Indeed, P ( 2 ) may intersect SN, or equivalently, 1224 C Budd and F Dun the set 9 may intersect SN", more than once. This is also indicated in figure 5.4.
The global effect of the curvature of W 3 on the subset P o means that, the velocity of the image of a point ( 4 , U) in 9 will increase with U if U is small, and will decrease with U for the larger values of U. As a consequence, the set PN+'@) lying in D will also be curved and the image of this set in 2 will retain this curvature. The effect on the map f(s) is for it to have a lower gadient or indeed a local minimum in the set Jw This is clearly evident in figure 3.8 showing the related map F(u) for U = -0.3375. The reduced gradient of f(s) in this subset means that it is now possible for a fixed point to be stable, and this is indeed what we observe. The windows of periodic behavior and the decrease in the number of impacts in these windows as E increases can both be explained in terms of values of E for which the curve Sk intersects the local maximum of W3, as the behaviour of the system must change qualitatively at such points. If we assume (very approximately indeed but consistently with the numerical observations!), that the maximum point of W 3 lies on a n extension of W obtained from the parameterization (4.2) and let cmax denote the value o f t at this point and if we assume further that the points of intersection of SK with Ware described (again approximately) by the formula 4.3, then the critical values of e = ek at which the behaviour changes are given by e&-'l ) / k so that ek+Jek = r if k is not too small. This is in agreement with the numerical observations. This roughly accounts for the numerical observations (i) and the observations (iii) can be explained in terms of the usual behaviour of continuous dynamical systems.