Dynamics of Multiple-Degree-of-Freedom Oscillators With Colliding Components

by examining the response of vibration absor bers with elastic stops, as well as the response of other non-linear oscillators with param eters leading to 2:1 and 3:1 internal resonance conditions. At the same time, the results illustrate some interesting features of the dynamics of these systems: namely, the absorber with stops is shown to possess superior technical characteristics. It provides a way in which to design non-linear absorbers and to suppress vibration levels over broader forcing frequency ranges than the conventional absorber, on both sides of the original resonant frequency. However, stability considerations in conjunction with numerical integration reveal that high amplitude beating and chaotic response may arise near the original reson ance, because the coexisting periodic response may become unstable due to Hopf bifurca tions. Next, similar phenomena are observed for other piecewise linear oscillators exhibiting 2:1 and 3:1 internal resonances. In these cases, zones of unstable periodic solutions are also observed in the vicinity of the middle of the primary resonances of the system. The resulting motions are again of much higher amplitude than that of the coexist· ing unstable n = 1 response. Several unstable zones may exist, generated by either Hopf or period doubling bifurcations. It is shown that the Hopf bifurcation of subharmonic solu tions results in multiple invariant circles. The coexistence of harmonic with large-period periodic and quasi-periodic responses is also presented. In all cases, similarities and differ ences in the response of systems possessing continuous non-linearities are noted.

by proceeding along the lines of the analysis for similar single-degree-of-freedom oscillators [13,14]. This type of analysis exploits the linearity of the equations of motion within each time interval where the system possesses constant characteristics. As expected [10], this extension is not a trivial task. A key feature that permits this generalization is the derivation of a suitable expression for the form of the solution of the linear problem.
The basic steps of the present analysis are briefly outlined in the following three sections. First, a suitable solution form is obtained for a general discrete linear system under harmonic excitation. This solution form then presents the basis for determining exact periodic motions of the non-linear system, as well as their stability. The analysis developed can be applied in a number of different areas. Here, the validity and effectiveness of the method is tested by applying it to two classes of important engineering problems. In this way, useful information is also obtained for the dynamics of the examples considered, which are vibration absorbers with elastic stops and other oscillators under internal resonance. The results of the numerical study reveal that some parameter combinations lead to a periodic or quasi-periodic response, while others lead to chaotically modulated motions. Comparison with the response of oscillators with continuous non-linearities is also made.

MECHANICAL MODEL AND SOLUTION FORM
The system examined is modeled as a general linear multiple-degree-of-freedom oscillator, including a component carrying elastic stops with damping and restoring forces, expressed by Fi . {ci+kx, (x, x)= c'x+k'x+(k-k')xs. x>xs.
In equation (1), x is the relative displacement of the oscillator carrying the stops with respect to the component of the system these stops collide with. An example of such a system with c =eN, c' =eN+ Cs, k = kN and k' = kN+ ks is presented in Figure 1. In general, the position of the stops within the overall system can be arbitrary.  where 4> is an unknown phase angle introduced for the purpose of facilitating the search of periodic solutions. Since the function F(x~ x) is piecewise linear, the exact form of the analytical solution of the system can be obtained for x ~ Xs and x < Xs. Before doing this, the displacements are normalized by a characteristic length, Xc, i.e., y 1 = xtf Xc. This length can be the gap length Xs, when Xs is non-zero. In addition, a normalized time is chosen according tor= rot. Then, the governing equations in the time intervals where the stiffness and damping remain constant can be cast in the linear form My+ Cy+ Ky=f('r) + k, (2) where M, C and K are the inertia, damping and stiffness matrices, respectively, while k is a vector including the constant terms of equation (1). Therefore, the general solution in the interval considered is y=yh+Yp· The homogeneous part, Yh• is expressed in terms of the eigenvalues s, = p 11 + iu 11 and the corresponding eigenvectors fn =ex,+ ip 11 , determined by where N is the total number of degrees of freedom. Note that the solution of the above eigenvalue problem is equivalent to that obtained by converting equation (2) to the classical first order 2N-space form. For oscillatory systems, the s,'s and y,'s appear in complex conjugate pairs. Making use of this, Yh is expressed in terms of real quantities as and combining all the above, the solution is finally expressed in the form where ~=(:).

PERIODIC STEADY STATE VIBRATIONS
Based on the solution form (4), a methodology is presented in this section, which leads to exact, n-periodic solutions of the system examined in a systematic and efficient way. These solutions are characterized by a single change of the damping and stiffness parameters within a response cycle. Since the system is non-linear, other vibrational modes are also possible. However, the motions considered are expected to be the most dominant. Then, the form of the solution sought is dictated by the form of u = y N, which is depicted in Figure 2 : namely, the solution is zl ( rl) in 0~ rl = r ~ rr and Zz( r z) in 0~ r 2 = r-rr ~ r!' where rr is the first contact time and r! = 2nn -rr . u Figure 2. The solution form of u = y N.
The determination of the desired solutions is based on the fact that they satisfy the following set of periodicity and matching conditions : The above represent a set of 4N + 2 transcendental equations for the 4N + 2 unknowns of the problem (i.e., y 1 , y 2 , t ~ and <P ). However, generalizing a solution procedure developed for similar but single-degree-of-freedom oscillators [ 13,14], the task of solving these equations can eventually be reduced to the solution of a single transcendental equation, as explained next.
First, application of equation (5) leads to (7) where w = [yl y 2 ]T. Note that for the solution form of Figure 2, ~1 =0. Next, for a given value of rr ' solution of equation (7) provides the vector w in the form w=wc cos q,+ws sin t/J+wo. (8) Then, application of equation (6) leads to two algebraic equations of the form bs sin q,+bc cos <P=bo. (9,10) Solving equations (9) and (10) for sin lfJ and cos t/J and applying a trigonometric identity, the correct value of -rf will give rise to (11) The strategy leading to a complete numerical solution of the problem is based on the observation that all the matrices and vectors of equation (7), as well as the a's and b's in equations (9) and (10), are functions of the system parameters and -rf only. Moreover, 0 < -rt < 2Hn. Then, the solution procedure starts by assigning values to -rt within the range from 0 to 21rn. This permits the determination of w from equation (7) and the evaluation of the corresponding value of G from equation ( 11 ), independently of the value of t/J. In doing so, the bisection method is applied in locating a possible solution for rr, when a sign change of G is detected between subsequent trial values of rf. Determination of rf allows the calculation of t/J from equations (9) and (10), the evaluation of the w from equation (8) and therefore the solution histories of Zt and z2 from equation (4). The above methodology can be modified and applied to other important cases, including systems with symmetric characteristics.

STABILITY OF PERIODIC RESPONSE
The information about the stability properties of a steady state motion is important because only stable solutions are realizable in practice. Moreover, the way in which a solution loses stability presents useful information in searching for other possible motions that may coexist. For piecewise linear systems, the classical methodologies are not applicable in performing the stability analysis of located periodic solutions. Here, a method employed before for other similar oscillators with one and two degrees of freedom will be generalized and applied.
Following reference [12], the essential information about the stability of a solution of a piecewise linear oscillator is obtained by just analyzing the propagation of arbitrary perturbations of the starting conditions within a single interval of two successive contact times. This actually implies that the steps of the analysis that follows are general and can be applied to solutions of piecewise linear systems of arbitrary form.
The stability procedure starts by letting .1 r 1 = .1 rt-.1 r3 , where .1 denotes a small perturbation of the variable it is applied on. Next, consider an exact periodic solution z 1 and a perturbed solution z 1 , obtained by perturbing the initial conditions of z 1 slightly (see Figure 2). Then, (12) where the vector Liz1 represents the deviation of Zt from Zt at time rt. Note that the corresponding variation of u1 is zero and, therefore, its position in .1z 1 is filled by the time variation .1 -rf instead. Application of equation (12) yields 4'.( rr + .1 r,)yt(Zo+ LiZo) + Jl( rr + .1 rt) cos 4> + v( rr + .1 tf} sin l/J = z.( rf) +Liz., (13) where the vector Li.zo involves the deviation from the correct initial conditions, Zo=z.(O), leading to z 1 • These initial conditions are determined by letting r 1 =0 in equation (4), once a periodic motion has been located by the analysis of the previous section. Next, expanding the left side of the last equation and keeping up to first order terms leads to tPt( rt)r.Lizo+ dJ.( rf)y.L1 '• +cit( rr)Li -rf =Liz., where q 1 ( r 1 ) = p., ( r 1 ) cos <P + v, ( r 1 ) sin <j), and the matrix r 1 is formed by adding the vector cl> ! 1 (0) q 1 (0) in the last column of the matrix cl> j" 1 (0). The zeroth order terms drop out because they satisfy equation (5). Then, the last equation can be put in the form with Q, = Rlt 1 R10. The matrix R 11 is formed by adding the vector -[ <i > 1 ( rf)y 1 +tit( rt)] in the last column of the identity matrix of order 2N x 2N, while R 10 is obtained by subtracting the vector <h1( rf)y, from the last column of the matrix cl>,( rt)r, .
Proceeding in a similar fashion in the interval 0 < r 2 < r~ yields where, again, the matrix Q 2 is found in terms of the system parameters only. Then, the relation between the perturbations at the beginning and at the end of the first response period is expressed by (14) where n = Q 2 Q 1 • Finally, after m response cycles, the perturbations from the correct initial conditions will be related to first order with the original perturbations by Therefore, if A. is the eigenvalue of ll with the largest modulus, the effect of the original error will vanish with time, provided that Ill < I. In such cases, the periodic solution examined is asymptotically stable. When IA.I > 1, the solution is unstable. Finally, when IA-1 = 1, bifurcations occur, resulting in qualitative changes of the system dynamics [15].
According to the results presented in the next section, three types of bifurcations were encountered in this study. The first type is a saddle-node bifurcation: it occurs with A.= 1 and results in the classical jump phenomena. The second type happens with ).. = -I and is associated with period doubling of the response. The third type is the so-called Hopf bifurcation, it occurs when a pair of complex eigenvalues has modulus one and results in an amplitude modulated response of the system. Note that if the damping coefficients are non-negative, the Hopf bifurcation is excluded for single degree of freedom oscillators [7,12], but it occurs for systems with N';?!::2 .

APPLICATIONS
This section verifies the analysis developed by presenting results for systems with two degrees of freedom (the system shown in Figure 1 with N = 2). The mass m 1 is subjected to a forcing F 1 cos r , while where The numerical results obtained for this system are presented in the following two subsections. In all cases considered, y 0 = -1. For other systems, the only difference in the analysis will be in the form of M, C, K, f and k.

VIBRATION ABSORBERS WITH ELASTIC STOPS
The idea of using stops in vibration absorption was presented in earlier studies [16,1]. The main interest there was to limit excessive amplitudes of vibration of the auxiliary mass of the absorber. Here, the analysis developed is applied in order to show some further enhancements and some limitations in the performance of the absorber with stops.
First, in Figure  originally in pretension and loses partial contact at y 0 and, as a consequence, presents softening characteristics, while the second non-linear oscillator possesses hardening characteristics. The effectiveness of the stops is obvious, when comparing the response diagram of the linear absorber (thin line) with the softening/ hardening absorber in the lower/ higher forcing frequency ranges, respectively. These benefits were expected from previous work on absorbers with continuous non-linearities (e.g., [17][18][19]). In fact, the present setup provides a practically feasible way of designing non-linear absorbers. However, in analogy to the non-linear absorber with continuous characteristics, the non-linearities also introduce some potential dangers, as explained next. In Figure 3, the unstable branch of the softening oscillator is generated with ).. = 1 bifurcations, associated with jump phenomena of the response at m = 0· 789 and 0· 806, and coexists with two other stable periodic solutions. Similar bifurcations are captured for the hardening system at m = I· 244 and 1· 259. In addition, the latter system exhibits period doubling (l= -1) bifurcations (at 0·805 and 0·816) as well as Hopf bifurcations (at 1·034 and 1·061). Response diagrams obtained by raising the non-linearity level are presented in Figure  4. This is done by increasing p' to 1· 732 and 4, which corresponds to stops with twice and 15 times the original stiffness of the absorber (approaching conditions of impact dampers), respectively, while maintaining 'i=0·02 and the values of the other parameters as before.
All the bifurcations of the lower resonance branches occur with ).. = -1, while ).. = I gives rise to the unstable motions of the higher resonance branches. Branch C is common to both diagrams. From the last plot, the most important observation is that increasing the non-linearity level causes a significant expansion of the zone of the unstable solutions originated by Hopf bifurcations near the original resonance ro = I. This is in accordance with results for absorbers with continuous non-linearities [ 18,19] and gives rise to rich but dangerous dynamics within that frequency range. For example, in Figure 5   (f) ciH ·0250.
As expected, a quasi·periodic response, represented by the invariant circle of Figure   5(a), appears at ro = 1·018, just after the Hopf bifurcation value of 1·017. This response becomes phase locked occasionally, as shown in Figure S(b) for the n = 13 solution obtained at 1·022. As the value of ro moves further away from the Hopfbifurcation value, the invariant circle expands, splits to n = 13 attractors ( Figure 5( c)) and distorts ( Figure  5(d)), according to a scenario explained in reference [20]. Then, following a single torus doubling ( Figure 5(e)), the system settles into chaotic response ( Figure 5(f)) without any further observable torus doubling (see references [21] and [22] for similar results and explanations).
To examine the effect of these motions on the performance of the absorber, in Figure   6(a) is shown the steady state response history of Y• at m = 1·018, during 95 forcing cycles. Clearly, the response involves beats, which is in agreement with previously reported experimental results [2,3]. In Figure 6 shows that caution should be exercised in choosing the absorber parameters properly. The present analysis can then be used as a guide in adjusting the parameters and eliminating these unstable responses.

OSCILLATORS WITH INTERNAL RESONANCE
The systems examined above have "equivalent" linearized counterparts with close natural frequencies. In this section, it is shown that some of the phenomena observed for the non-linear absorber can also occur in other oscillatory systems, when their parameters are chosen in ways that lead to internal resonance conditions [23][24][25][26]. A number of some new and interesting features also appear and are analyzed.
First, in Figure 7 is presented the response diagram for a system with parameters p = fJ = 0· 5, p = 0·9, p' = 0·65, '• = 0·01 and ' 2 = ' 2 = 0·02. This combination of parameters leads to 2: 1 internal resonance conditions. First, note that there exist two resonance peaks around ro = 0· 65 with very close frequencies, while the peak of the resonance branch around 1· 3 consists of unstable solutions, generated by A.= -1 bifurcations at 1· 203 and 1· 348. As expected, these bifurcations give rise to a periodic response with a period equal two times the forcing period. To illustrate this, n = 2 solutions arising in that frequency range are also determined and shown. In addition, there is a branch of n = 1 unstable solutions in the range from ro =0·963 to 0·990, originated by a Hopf bifurcation. Again, this branch is located around the middle of the primary resonance frequencies. However, there is a qualitative difference between the response of the present system and that of the oscillators analyzed in the previous section: namely, varying the last set of parameters may result in unstable ranges near the middle of the primary resonances, generated by period doubling and not by Hopf bifurcations. For other sets of parameters, several unstable intervals may appear in the same range, some generated by period doubling and the others generated by Hopf bifurcations. For example, decreasing the damping ratios of the last example to ' 1 = ' 2 = ' 2 = 0· 001, apart from an order of magnitude change of the response amplitude and a slight bilateral expansion of the Hopf generated unstable solutions (which now occupy the range of w from 0·955 to 0·998), an additional unstable range of n = 1 solutions is generated with A. = -1 between 0·828 and 0·865.
Similar observations are made for the response of the following example, with parameters JJ = P = 1, p=0·4, p' =0· 35 and '• =' 2 ='2=0·001, leading to 3:1 internal resonance conditions. The response diagram for this system is shown in Figure 8. Again, an unstable range of n = 1 solutions appears through A.= -1 bifurcations at ro = 0· 685 and 0· 732. The n = 2 solutions existing in that frequency range are also presented in order to illustrate that  their amplitude is larger (by more than an order of magnitude near resonance) than the amplitude of the coexisting n = 1 motions. Moreover, the n = 2 solutions undergo saddlenode bifurcations, giving rise to a very narrow frequency range in which two stable and one unstable n = 2 solutions coexist. Finally. note the multiple peaks for both then= 1 and n = 2 branches. As for the 2: 1 case, changing the parameters may result in multiple frequency ranges involving unstable response generated by either A.= -1 or Hopf bifurcations. For example, in Figure 9 is shown the response diagram obtained for the last example between m = 0·8 and 1·2 by increasing the value of p to 0·9. Both n = 1 (index 1) and n = 2 (index 2) solutions are shown. The unstable branches denoted by H and D are generated by Hopf and period doubling bifurcations, respectively. In fact, in some cases the bifurcation takes place with a pair of complex conjugate eigenvalues equal to one, while a real eigenvalue has a value very close to -1, indicating conditions close to codimension two bifurcations [15].
To present a small but representative sample of some of the interesting dynamics encountered in that range, in Figure 10 is shown a sequence of Poincare plots (y 2 ,j 2 ) obtained by direct integration for the parameters of the last example. For that system, a Hopf bifurcation of the n = 2 solution occurs at m = 1· 1292. Poincare plots at m = 1· 130 and 1· 129 are shown in Figure 1 0( a). Obviously, the n = 2 solution, represented by the two heavy dots, is replaced by two invariant circles, following the bifurcation. Moving to lower values of m, the solution finally locks to ann= 8 solution, as shown in Figure lO is also determined at the same frequency. The interesting feature of this solution is that it consists of nine invariant circles, possibly originated by a Hopf bifurcation of an n = 9 solution. When m is varied, each of these invariant circles undergoes changes similar to those of the single circle of Figure 5(a). The behavior encountered in this section presents similarities and differences with the behavior of systems with asymmetric but continuous non-linearities [24]. Here, in addition to the unstable zones generated by the Hopf bifurcation between the primary resonant frequencies, there may exist zones that originated with A,= -1 bifurcations. This last type of bifurcation does not have an analog in the analysis of the steady state solutions captured in reference [24]. The reason is that the stability of those solutions is governed by linearized differential equations (flows), while here it is governed by difference equations (maps).

SYNOPSIS AND CONCLUSIONS
An investigation is carried out on the long-time response of a class of multiple-degreeof-freedom oscillators with strongly non-linear characteristics. The non-linearity of the system is modeled as a component with piecewise linear damping and stiffness properties. The solution procedure exploits the piecewise linear nature of the governing equations. Suitable solution forms are first developed for the response within the intervals where the system characteristics remain constant. Then, exact, n-periodic, single crossing steady state motions are determined by imposing a set of periodicity conditions. For an N-degree-offreedom system, this results in a system of 4N + 2 transcendental equations. However, generalizing previous work on similar single-degree-of-freedom systems, the numerical effort is reduced to finding the solution of sets of linear equations and eventually to a single transcendental equation for the crossing time. Finally, the analytical procedure is complemented by a stability analysis which is appropriate for solutions of systems with arbitrary piecewise linear characteristics. As in the case of single-degree-of-freedom oscillators, this task is reduced to examining the propagation of errors in starting conditions within a single linear motion interval.
The present analysis provides efficient means for performing parametric studies in a large number of systems encountered in engineering practice. Here, the applicability and effectiveness of this analysis is demonstrated by examining the response of vibration absorbers with elastic stops, as well as the response of other non-linear oscillators with parameters leading to 2:1 and 3:1 internal resonance conditions. At the same time, the results illustrate some interesting features of the dynamics of these systems: namely, the absorber with stops is shown to possess superior technical characteristics. It provides a way in which to design non-linear absorbers and to suppress vibration levels over broader forcing frequency ranges than the conventional absorber, on both sides of the original resonant frequency. However, stability considerations in conjunction with numerical integration reveal that high amplitude beating and chaotic response may arise near the original resonance, because the coexisting periodic response may become unstable due to Hopf bifurcations. Next, similar phenomena are observed for other piecewise linear oscillators exhibiting 2:1 and 3:1 internal resonances. In these cases, zones of unstable periodic solutions are also observed in the vicinity of the middle of the primary resonances of the system. The resulting motions are again of much higher amplitude than that of the coexist· ing unstable n = 1 response. Several unstable zones may exist, generated by either Hopf or period doubling bifurcations. It is shown that the Hopf bifurcation of subharmonic solutions results in multiple invariant circles. The coexistence of harmonic with large-period periodic and quasi-periodic responses is also presented. In all cases, similarities and differences in the response of systems possessing continuous non-linearities are noted.