# Stable periodic vibroimpacts of an oscillator

N. Popplewell, Charusheel Barat, K N Mclachlan

► To cite this version:

N. Popplewell, Charusheel Barat, K N Mclachlan. Stable periodic vibroimpacts of an oscillator. Journal of Sound and Vibration, 1983, 87, pp.41 - 59. 10.1016/0022-460X(83)90438-8. hal-01510833

HAL Id: hal-01510833

https://hal.science/hal-01510833

Submitted on 20 Apr 2017

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
STABLE PERIODIC VIBROIMPACTS OF AN OSCILLATOR

N. POPPLEWELL AND C. N. BAPAT

*Department of Mechanical Engineering*

AND

K. MCLACHLAN

*Department of Civil Engineering, University of Manitoba, Winnipeg, Canada R3T 2N2*

The asymptotically stable vibrations of a loaded oscillator colliding periodically with a rigid mass are described. Comparison of the numerical results with the few existing examples is encouraging but inconclusive. Better overall agreement is demonstrated with fairly comprehensive measurements from a specially built experimental rig. Impact motions are shown to be very sensitive to small fluctuations in the clearance between masses and the stiffness and loading of the oscillator near its linear, or collisionless, resonant frequency. The rigid mass is quite an effective damper at or just above this frequency condition.

# 1. INTRODUCTION

A general theory to describe the stable periodic motions of two colliding rigid bodies called an impact-pair has been developed and verified experimentally in a companion paper [1]. The theory will be extended by letting one of the bodies have flexibility and damping in order to assess resonant effects. Such effects are likely important in high speed mechanisms or machinery where clearances arise from tolerances, wear, etc. On the other hand, an extra impacting body may be introduced deliberately into a resonant system to reduce excessive vibration excursions. The device is normally termed an "impact damper" and can be quite effective despite the high transient accelerations generated by the impacts [2-21]. The basic theory of these applications will be checked against previous, rather restricted, theoretical results and new, exhaustive experimental data.

A closed form solution has been derived for impacts consistently happening twice for every cycle of a sinusoidal load [2-5, 7, 8, 10, 11, 14, 22, 23]. The impacts have been assumed generally instantaneous and distributed primarily evenly in time. Masri [7, 13] and the Kobrinskii's [23] subsequently applied the concept of error propagation in difference equations to ascertain whether these periodic motions were asymptotically stable or not. These approaches are fundamental but the practical application of the solutions is quite restricted. For example Sadek [6] and Dittrich [9] suggested from limited experimental and theoretical studies that uneven temporal distributions are most likely at frequencies close to where the collisionless, flexible system is resonant. Sadek also inferred after careful experimentation that "equally spaced impacts hardly ever occur" even away from this resonance. Therefore it is not surprising that Sadek [6] and Masri [13] attempted an extension to unevenly spaced impacts. Two impacts per loading
cycle were emphasized; only a single, probably analog computed illustration of three impacts was given by Masri alone. The lone three impact example may well be the consequence of a slowly converging and, hence, computationally time consuming iteration scheme for non-linear equations. An efficient algorithm will be used here for a comparable set of linearized equations to facilitate the computation of stable periodic motions involving any number of impacts. Numerical examples however will be confined to five or fewer impacts in each loading cycle. The theory basically is an elaboration of Masri's work which is itself based on the concept of a coefficient of restitution and the conservation of linear momentum. Previous experiments suggest that answers may not be unique in those instances when the theoretically neglected transients cannot be ignored [1].

## 2. GENERAL THEORY

The idealized vibroimpact device is shown in Figure 1. Its primary system consists of a linear spring with stiffness $M$, a viscous dashpot having damping constant $C$, and mass

<figure><img src="image_2.png" /><figcaption>Figure 1. Single degree of freedom system with a second impacting mass.</figcaption></figure>

$M$ excited by the external harmonic load, $F_0 \sin \Omega t$. The secondary system is composed of a rigid mass, $m$, which can move uniaxially in a slot inside the usually much heavier mass $M$. The supposed frictionless motion of $m$ is instigated by collisions with $M$ which occur intermittently because of clearance $d$. Any one impact, typically $i$, is assumed, reasonably for metals [1, 22, 23], to be instantaneous and to be described by a coefficient of restitution, $R_i$. This coefficient may vary periodically in the general formulation but constant values, $R_1$ and $R_2$ corresponding to collisions on the right- and left-hand sides of the slot in Figure 1, will be stipulated in illustrative examples.

Equations of motion and their solutions will be developed in a manner similar to that in reference [1]. The main difference is that the conservation of linear momentum now has to be employed along with the definition of $R_i$ to accommodate the discontinuous velocities at a collision. Displacements, conversely, are invariably continuous and the sinusoidal load sustains motions of the primary and secondary systems which bear similar periodicity conditions to before. However, the periodicity is related more conveniently to the temporal behaviour of the load rather than the displacement of the primary mass. This relationship is illustrated in Figure 2. A motion's periodicity is shown in Figure 2(b) to be $T_0$ during which time the load in the particular example of Figure 2(a) has undergone two cycles. Meanwhile $N$ impacts between $m$ and $M$, represented by the spikes in Figure 2(c), have happened. Consequently this motion will be labelled $N$ impacts/two cycles (of the load) with the generalization being obvious. The motion is periodic so that spike $Y_{(N+1)}$ at instant $t_{(N+1)}$, equal to $(t_1 + T_0)$, is identical to $Y_1$, $Y_{(N+2)}$ at $t_{(N+2)}$, or $(t_2 + T_0)$, is identical to $Y_2$, etc. The traditional approach will be taken of considering the periodicity, $T_0$, the number of impacts, $N$, and their occurrence at instants $t_1, t_2, \dots, t_N$ as known. The load's amplitude and phase with respect to the first impact on the other hand will be treated together with the absolute displacement of $M$ at each of the $N$ collisions as unknown. Relationships between these $(N+2)$ unknowns will be developed next.
<figure><img src="image_3.png" /><figcaption>Figure 2. General periodic motion of N impacts happening during period $T_0$. (a) Load $F_0 \sin \Omega t$; (b) distribution in time of the N impacts; (c) presumed sequence of impacts.</figcaption></figure>

The equation of motion of mass M between impacts is

$$M\ddot{X} + C\dot{X} + KX = F_0 \sin \Omega t, \quad (1)$$

where $X(t)$ represents the absolute displacement of the primary mass M. A dot superscript indicates differentiation with respect to time. The solution of equation (1) can be obtained by superimposing free and forced motions [5, 24]. A typical solution between two consecutive collisions, say the $i$ and the $(i+1)$ impacts happening at $t_i$ and $t_{(i+1)}$ respectively, is

$$X(t) = \exp\left\{-\frac{\xi}{r}[\Omega t - \alpha_i]\right\}\left\{a_i \sin\left(\frac{\eta}{r}[\Omega t - \alpha_i]\right) + b_i \cos\left(\frac{\eta}{r}[\Omega t - \alpha_i]\right)\right\} + A \sin[\Omega t + \tau], \\ t_{ia} \le t \le t_{(i+1)b}, \quad i = 1, 2, 3, \dots, N. \quad (2)$$

The temporal origin is taken at $t_1$ so that $t$ in equation (2) is the convenient abbreviation of $(t-t_1)$. Quantities just after and just before an impact are represented by subscripts a and b respectively. Consequently $t_{ia}$, for example, is the instant immediately after the $i$th impact. Variables $\xi, \omega, r, \eta, \tau$ and A are given conventionally by

$$\xi = C/2(KM)^{1/2}, \quad \omega = (K/M)^{1/2}, \quad r = \Omega/\omega, \quad \eta = (1 - \xi^2)^{1/2}, \quad (3)$$

$$\tau = \Omega t_1 - \tan^{-1}\{2\xi r/(1-r^2)\}, \quad A = (F_0/K)\{(1-r^2)^2 + (2\xi r)^2\}^{-1/2}. \quad (4)$$

Remaining variables are defined as

$$a_i = (1/\eta)\{(1/\omega)\dot{X}_{ia} - \text{Arcos}[\alpha_i + \tau] + \xi b_i\}, \quad i = 1, 2, \dots, N, \quad (5)$$

$$b_i = X_i - A \sin[\alpha_i + \tau], \quad i = 1, 2, \dots, N, \quad (6)$$

$$X_i \equiv X(t_i), \quad \alpha_i = \begin{cases} 0, & i=1 \\ \Omega t_i, & i=2, 3, \dots, N \end{cases}. \quad (7)$$

$(\alpha_1$ is zero because the temporal origin has been taken as the instant of the first impact.) The displacement of M is continuous at collisions so that equation (2) remains valid then. Consequently, the unknown displacement of M at the N collisions within $T_0$, $X_i$, may be obtained by merely substituting their instants of occurrence, $t_i$ where $i = 1, 2, \dots, N$, into equation (2). For example, the displacement at the $(i+1)$th impact is

$$X(t_{i+1}) \equiv X_{(i+1)} = C_{1i}a_i + C_{2i}b_i + A \sin[\alpha_{(i+1)} + \tau], \quad i = 1, 2, \dots, N, \quad (8)$$
where, by employing equation (7),

$$C_{1i} = \exp \left\{ -\frac{\xi}{r} (\alpha_{i+1} - \alpha_i) \right\} \sin \left\{ \frac{\eta}{r} (\alpha_{i+1} - \alpha_i) \right\}, \quad (9)$$

$$C_{2i} = \exp \left\{ -\frac{\xi}{r} (\alpha_{i+1} - \alpha_i) \right\} \cos \left\{ \frac{\eta}{r} (\alpha_{i+1} - \alpha_i) \right\}, \quad i = 1, 2, \dots, N. \quad (10)$$

The $X_{(i+1)}$ and, from equation (4), A and $\tau$ are certainly not known in these last three expressions. The $\alpha_i$, $i = 1, 2, \dots, N$, may be calculated from the given driving frequency, $\Omega$, and impact instants, $t_i$, by using equation (7). Then the $C_{1i}$ and $C_{2i}$ may be evaluated by substituting these $\alpha_i$ and the determinable parameters of equation (3) into equations (9) and (10). Similar substitutions into equations (5) and (6) leave the final two variables in equation (8), $a_i$ and $b_i$, requiring $\dot{X}_{ia}$ other than the $X_i$, A and $\tau$. The discontinuity in the velocity of M at the $i$th impact may be obtained from the conservation of linear momentum and the definition of $R_i$.

Linear momentum conservation gives [5, 13, 23]

$$M\dot{X}_{ib} + m\dot{Z}_{ib} = M\dot{X}_{ia} + m\dot{Z}_{ia}, \quad i = 1, 2, \dots, N, \quad (11)$$

for the $i$th impact whereas the definition of the coefficient of restitution is

$$R_i = -(\dot{X}_{ia} - \dot{Z}_{ia})/(\dot{X}_{ib} - \dot{Z}_{ib}), \quad i = 1, 2, \dots, N. \quad (12)$$

The $\dot{Z}_{ib}$ and $\dot{Z}_{ia}$ are the absolute velocities of the secondary mass just before and after the $i$th impact, respectively. Equations (11) and (12) may be solved straightforwardly to give

$$\dot{X}_{ia} = \{(1 + \mu)R_i/(1 + R_i)\}\dot{Z}_{ib} + \{(1 - \mu R_i)/(1 + R_i)\}\dot{Z}_{ia}, \quad i = 1, 2, \dots, N, \quad (13)$$

$$\mu = m/M \quad (14)$$

being the (known) mass ratio of the secondary and primary systems. The velocity of $m$ must remain constant between two consecutive impacts because no forces are assumed to act upon it then. Consequently, its velocity $\dot{Z}_{ia}$ between the $i$ and $(i+1)$ impacts, for example, is simply the ratio of the absolute distance travelled by $m$ and the time elapsed between these two impacts. Mathematically,

$$\dot{Z}_{ia} = (Z_{(i+1)} - Z_i)/(t_{(i+1)} - t_i), \quad i = 1, 2, \dots, N. \quad (15)$$

Substituting equation (7) into equation (15) produces

$$\dot{Z}_{ia} = \frac{\Omega[X_{(i+1)} + Y_{(i+1)} - X_i - Y_i]}{\alpha_{(i+1)} - \alpha_i}, \quad Y_i = Z_i - X_i, \quad i = 1, 2, \dots, N \quad (16, 17)$$

being the relative displacement between $m$ and $M$. The $Y_i$ is introduced because impacts occur only when the relative displacement is $+d/2$ for a collision on the right side of $M$ in Figure 1 and $-d/2$ for a collision on the left. Hence, at impact,

$$Y_i = \begin{cases} d/2, & \text{collision on right of } M \\ -d/2, & \text{collision on left of } M \end{cases}, \quad i = 1, 2, \dots, N. \quad (18)$$

Consequently, the only unknowns in expression (16) for $\dot{Z}_{ia}$ are the $X_i$. Similarly, the velocity of the secondary mass just before the $i$th impact, $\dot{Z}_{ib}$, equals its velocity between the $(i-1)$th and $i$th impacts. Therefore,

$$\dot{Z}_{ib} = \frac{\Omega[X_i + Y_i - X_{(i-1)} - Y_{(i-1)}]}{\alpha_i - \alpha_{(i-1)}}, \quad i = 1, 2, \dots, N. \quad (19)$$
Equations (13), (16) and (19) may be combined to give

$$ \dot{X}_{ia} = \frac{(1+\mu)R_i\Omega}{(1+R_i)(\alpha_i - \alpha_{(i-1)})} [X_i - X_{(i-1)} + Y_i - Y_{(i-1)}] \\ + \frac{(1-\mu R_i)\Omega}{(1+R_i)(\alpha_i \upharpoonright_{13} - \alpha_i)} [X_{(i+1)} - X_i + Y_{(i+1)} - Y_i], \quad i=1, 2, \dots, N. \quad (20) $$

This last relationship may be substituted into equation (5) to produce an expression for $a_i$ in terms of the unknown $X_i$, $A$ and $\tau$. Consequently,

$$ a_i = \frac{(1+\mu)R_i\Omega}{(1+R_i)\eta\omega(\alpha_i - \alpha_{(i-1)})} [X_i - X_{(i-1)} + Y_i - Y_{(i-1)}] \\ + \frac{(1-\mu R_i)\Omega}{(1+R_i)\eta\omega(\alpha_{(i+1)} - \alpha_i)} [X_{(i+1)} - X_i + Y_{(i+1)} - Y_i] - \frac{Ar}{\eta} \cos[\alpha_i + \tau] + \frac{\xi}{\eta} b_i, \quad i=1, 2, \dots, N, \quad (21) $$

where, the $b_i$ of equation (6) are already in the desired form,

$$ b_i = X_i - A \sin[\alpha_i + \tau], \quad i = 1, 2, \dots, N. \quad (22) $$

Equations (9), (10), (17), (21) and (22) can be substituted into equation (8) to give $N$ linear, simultaneous equations.

A perusal of relation (21) should indicate that terms like $\alpha_{(i-1)}$ and $\alpha_{(i+1)}$ have not yet been defined at the two extreme values of $i$. They have to be obtained from the requirement for periodic motion. It has been seen already that periodicity implies $t_{(N+1)} = t_1 + T_0$, $t_{(N+2)} = t_2 + T_0$, ..., or, in general,

$$ t_{(N+i)} = t_i + T_0, \quad Y_{(N+i)} = Y_i, \quad X_{(N+i)} = X_i, \quad R_{(N+i)} = R_i, \quad i = 1, 2, \dots, N. \quad (23, 24) $$

Multiplying equation (23) by $\Omega$ and using equation (7) leads to

$$ \alpha_{(N+i)} = \alpha_i + \Omega T_0 = \alpha_i + 2\pi k, \quad i = 1, 2, \dots, N \quad (25) $$

if the load undergoes $k$ cycles in time $T_0$. Consequently, when $i$ equals $N$, $\alpha_{(i+1)}$ is simply $(\alpha_1 + \Omega T_0)$, or $\Omega T_0$ in view of equation (7). When $i$ is unity on the other hand, terms like $(\alpha_i - \alpha_{(i-1)})$ in relation (21) may be written as $(\alpha_{(N+i)} - \alpha_{(N+i-1)})$ by employing equation (25). These last expressions can be evaluated by using the known driving frequency, $\Omega$, the assumed periodicity, $T_0$, and presumed contact instants $t_i$, $i = 1, 2, \dots, N$. Terms such as $Y_{(N+i-1)}$ can be computed by further presuming the side of $M$ at which any particular impact happens. Then equation (18) can be used straightforwardly. All other variables except for the $X_i$, $i = 1, 2, \dots, N$, $A$ and $\tau$ are understood to be given or calculable from equations (3), (14) and (24). Consequently ($N+2$) unknowns exist. If trigonometric terms like $\sin[\alpha_{(i+1)} + \tau]$ in equations (8) and (22) are expanded and the resulting coefficients of the unknowns unified, the $N$ simultaneous equations can be shown to be [25]

$$ W_{1i}X_i + W_{2i}X_{(i+1)} + W_{3i}X_{(N+i-1)} + W_{4i}A \cos \tau + W_{5i}A \sin \tau = W_{6i}, \quad i = 1, 2, \dots, N. \quad (26) $$

The often cumbersome but calculable coefficients $W_{ri}$, $r = 1, 2, \dots, 6$ and $i = 1, 2, \dots, N$, are given more conveniently in the Appendix. The $X_i$ must still satisfy periodicity condition (24) so that only $X_1, X_2, \dots, X_N$ are actually involved in equation (26).

An additional $N$ linear simultaneous equations may be developed by considering the velocity rather than primarily the displacement of $M$. The velocity of $M$ just before the
(i + 1)th impact, for example, is simply the derivative of equation (2) at $t_{(i+1)b}$. An alternative expression may be found like that for the velocity just after an impact from linear momentum and the definition of the coefficient of restitution. Relating the two expressions produces the desired equations. From equation (2)

$$ 
\begin{aligned}
\dot{X}_{(i+1)b} = & \exp \left\{ -\frac{\xi}{r}(\alpha_{(i+1)} - \alpha_i) \right\} \left[ \left\{ a_i \cos \frac{\eta}{r}(\alpha_{(i+1)} - \alpha_i) - b_i \sin \frac{\eta}{r}(\alpha_{(i+1)} - \alpha_i) \right\} \frac{\eta\Omega}{r} \right. \\
& \left. - \frac{\xi}{r}\Omega \left\{ a_i \sin \frac{\eta}{r}(\alpha_{(i+1)} - \alpha_i) + b_i \cos \frac{\eta}{r}(\alpha_{(i+1)} - \alpha_i) \right\} \right] + A\Omega \cos[\alpha_{(i+1)} + \tau],
\end{aligned}
\quad i = 1, 2, \dots, N, \quad (27) $$

with equation (7) and $\dot{X}_{(i+1)b}$ rather than $\dot{X}(t_{(i+1)b})$ being employed. Solving the momentum and restitution equations (11) and (12) for $\dot{X}_{ib}$ rather than $\dot{X}_{ia}$ gives

$$ \dot{X}_{ib} = \{(R_i - \mu)/(1 + R_i)\}\dot{Z}_{ib} + \{(1 + \mu)/(1 + R_i)\}\dot{Z}_{ia}, \quad (28) $$

where $\mu$ is defined by equation (14). Expressions (16) and (19) for $\dot{Z}_{ia}$ and $\dot{Z}_{ib}$ may be substituted into equation (28) to yield

$$ 
\begin{aligned}
\dot{X}_{ib} = & \left(\frac{R_i - \mu}{1 + R_i}\right) \frac{\Omega}{(\alpha_i - \alpha_{(i-1)})} [X_i - X_{(i-1)} + Y_i - Y_{(i-1)}] \\
& + \left(\frac{1 + \mu}{1 + R_i}\right) \frac{\Omega}{(\alpha_{(i+1)} - \alpha_i)} [X_{(i+1)} - X_i + Y_{(i+1)} - Y_i], \quad i = 1, 2, \dots, N.
\end{aligned}
\quad (29) $$

Now $\dot{X}_{(i+1)b}$ is needed rather than $\dot{X}_{ib}$. It is obtained by replacing every $i$ in equation (29) with $(i+1)$. Identifying the resulting relationship to equation (27) leads to the following simultaneous equations [25]:

$$ V_{1i}X_i + V_{2i}X_{(i+1)} + V_{3i}X_{(i+2)} + V_{4i}X_{(N+i-1)} + V_{5i}A \cos \tau + V_{6i}A \sin \tau = V_{7i}, \quad i = 1, 2, \dots, N, \quad (30) $$

after expanding trigonometric terms and accounting for periodicity in the manner needed to generate equation (26). The computable coefficients $V_{ri}$, $r = 1, 2, \dots, 7$ and $i = 1, 2, \dots, N$, are again detailed for a given $i$ in the Appendix. The $(N+2)$ unknowns of $X_i$, $i = 1, 2, \dots, N$, $A$ and $\tau$ are common to both equations (26) and (30) because of the periodicity condition (24). Consequently $2N$ linear simultaneous equations have been formed in $(N+2)$ unknowns$^\dagger$.

Exact closed form solutions are possible when the number of equations and unknowns are identical ($N$ equal to two). Then the $X_i$, $i = 1, 2, \dots, N$, can be eliminated completely leaving only $A$ and $\tau$ as unknowns [25]. The two resulting equations are independent except in the special, essentially symmetric case when coefficients of restitution are identical and the durations between any three consecutive impacts are the same. This situation is analogous to that observed previously for exclusively rigid impacting mechanisms [1]. Consequently the same procedure as before was adopted: i.e., assume some value for $A$ then calculate the ensuing $\tau$. More generally however, the IMSL subroutine LEQT1F was used to solve the four independent equations (26) and (30) numerically by utilizing Gaussian elimination [26]. This course was felt more expedient than solving the equations analytically to produce extremely lengthy expressions. A numerical tech-

<footnote>$^\dagger$Alternatively, the $A$ and $\tau$ may have been assumed rather than the $\alpha_i$ and the $2N$ unknowns of $X_i$ and $\alpha_i$, $i = 1, 2, \dots, N$, calculated, maybe somewhat more precisely, from the $2N$ equations. However, the resulting equations are linear no longer so that numerical computations are very likely to be much more time consuming.</footnote>
nique appears virtually essential on the other hand when N exceeds two. Then the least squares subroutine LLSQF of the IMSL library [26] was employed. Values of the unknowns were finalized when the square of the difference between the left and right sides of relations (26) and (30), $\epsilon^2$, was reduced so that $\epsilon$ was less than 3% of these values. (The 3% was decreased to 1.5% with virtually no effect in several arbitrarily chosen examples.) All computations were performed throughout the present work with double precision arithmetic on an AMDAHL V7 digital computer.

Periodic solutions are grouped conventionally into asymptotically stable or unstable subsets. A periodic motion is asymptotically stable if, after a small perturbation, the original motion is restored eventually [7, 13, 23]. Consequently asymptotic stability would enable machinery performance to be predictable [22, 27-29]. On the other hand, asymptotic instability implies that a small perturbation will change some aspect of the periodic motion for all time. The asymptotic stability of periodic motions was determined by the Kobrinskii's [23] and Masri [7, 13] for contact surfaces having a single constant coefficient of restitution. Small perturbations at or just after the first impact were related to corresponding changes at or just after the second, third, etc., impacts. This technique was adapted in a companion paper to a rigid impact-pair with different, maybe periodically varying, coefficients of restitution. A similar, somewhat minor extension is possible here. The absolute displacement and velocity of M are no longer prescribed here so they need perturbing along with the previously changeable Z and $\dot{Z}$. Some of the presumed collision instants, and hence the $\alpha_i$, may be modified in accordance with the equations of motion and periodicity conditions. However, these adjustments must be transitional for asymptotically stable motion because the timing sequence must revert eventually to the claimed periodic form. Masri [13] adopted a slightly different approach but with the same end results. He replaced Z, which equation (16) indicates depends upon $\alpha_i$, by the $\alpha_i$ of equation (7). Then the initial perturbation and subsequent changes, each denoted by $\Delta$, were shown to be related by

$$ \begin{Bmatrix} \Delta X_{(N+1)a} \\ \Delta \dot{X}_{(N+1)a} \\ \Delta \dot{Z}_{(N+1)a} \\ \Delta \alpha_{(N+1)a} \end{Bmatrix} = [P] \begin{Bmatrix} \Delta X_{1a} \\ \Delta \dot{X}_{1a} \\ \Delta \dot{Z}_{1a} \\ \Delta \alpha_{1a} \end{Bmatrix}, \quad [P] = [P_{Na}][P_{(N-1)a}] \cdots [P_{1a}] \cdots [P_{1a}], \quad (31, 32) $$

if terms involving a single perturbation are significantly larger than expressions containing products of perturbations. Nomenclatures of previous references have been adjusted to be compatible with earlier formulations in this section. The $a$ in a subscript for instance still implies an instant immediately after an impact. Different coefficients of restitution have been accommodated in the elements of the typical component matrix $[P_{ia}]$ given by equations (A5) through (A7) of the Appendix. These coefficients and variables like $X_i$ and $Y_i$ must fulfil periodicity condition (24) and the $\alpha_i$ must similarly comply with equation (25). In addition, the $Y_i$ must satisfy equation (18). All elements of the component matrices may be calculated once the unknown $X_i$, $i = 1, 2, \dots, N$, $A$ and $\tau$ have been computed in the fashion detailed previously. Then $[P]$ may be formed according to equation (32). If all the moduli of the eigenvalues of $[P]$ are strictly less than unity, then a periodic motion is asymptotically stable [13]. Eigenvalues were computed by employing Hessenberg's method in the readily available EIGRF of the IMSL library [26].

## 2.1. COMPARISON WITH PREVIOUS THEORETICAL RESULTS

Computations were performed by using the theory of the last section in order to compare the results with more limited theoretical data available for several particular
<figure><img src="image_4.png" /><figcaption>Figure 3. Comparison of predicted two equispaced, two and three unequispaced impacts/cycle stability zones with the previously most comprehensive results at $d=37(F_0/K)$ alone. $\dagger$, Imples equispaced, e, or unequispaced, u, impacts/cycle. $\xi=0.005$; $\mu=0.042$; $R_1=R_2=0.75$.</figcaption></figure>

<caption>TABLE 1</caption>

A numerical comparison of previous [13] and present theoretical results for two equispaced and two unequispaced impacts/cycle of the external load

Parameters of primary and secondary system

$$ \xi = 0.005, r = 0.94, \mu = 0.042, R_1 = R_2 = 0.75, d/(F_0/K) = 37 $$

<table><thead><tr><th rowspan="2"></th><th colspan="2">2 Equispaced†</th><th colspan="2">2 Equispaced</th><th colspan="2">2 Unequispaced</th></tr><tr><th>Present results</th><th>Previous results</th><th>Present results</th><th>Previous results</th><th>Present results</th><th>Previous results</th></tr></thead><tbody><tr><td>α<sub>2</sub>/π</td><td>1</td><td>1</td><td>1</td><td>1</td><td>1.37</td><td>1.37</td></tr><tr><td>Ωt<sub>1</sub>/π</td><td>0.68</td><td>0.68</td><td>-0.53</td><td>-0.53</td><td>0.52</td><td>0.52</td></tr><tr><td>τ/π</td><td>0.66</td><td>0.65</td><td>-0.56</td><td>-0.55</td><td>0.50</td><td>0.49</td></tr><tr><td>X<sub>1</sub></td><td>16.17</td><td>16.21</td><td>-5.10</td><td>-5.07</td><td>13.65</td><td>13.76</td></tr><tr><td>X<sub>2</sub></td><td>-16.17</td><td>-16.21</td><td>5.10</td><td>5.07</td><td>-9.31</td><td>-9.37</td></tr><tr><td>Z<sub>1</sub></td><td>34.67</td><td>34.71</td><td>13.41</td><td>13.42</td><td>32.15</td><td>32.26</td></tr><tr><td>Z<sub>2</sub></td><td>-34.67</td><td>-34.71</td><td>-13.41</td><td>-13.42</td><td>-27.81</td><td>-27.87</td></tr><tr><td>X&#770;<sub>1a</sub></td><td>-2.22</td><td>-2.22</td><td>-0.86</td><td>-0.85</td><td>5.47</td><td>5.56</td></tr><tr><td>X&#770;<sub>2a</sub></td><td>2.22</td><td>2.22</td><td>0.86</td><td>0.85</td><td>9.81</td><td>10.03</td></tr><tr><td>X&#770;<sub>1b</sub></td><td>-3.96</td><td>-3.96</td><td>-1.53</td><td>-1.53</td><td>3.73</td><td>3.80</td></tr><tr><td>X&#770;<sub>2b</sub></td><td>3.96</td><td>3.96</td><td>1.53</td><td>1.53</td><td>11.66</td><td>11.78</td></tr><tr><td>Z&#770;<sub>1a</sub></td><td>-20.74</td><td>-20.77</td><td>-8.02</td><td>-8.03</td><td>-13.09</td><td>-13.10</td></tr><tr><td>Z&#770;<sub>2a</sub></td><td>20.74</td><td>20.77</td><td>8.02</td><td>8.03</td><td>28.48</td><td>28.69</td></tr><tr><td>X<sub>max</sub>‡</td><td>-16.67, 16.67</td><td>-16.72, 16.72</td><td>-5.18, 5.18</td><td>-5.16, 5.16</td><td>-15.17, 14.77</td><td>-15.32, 13.79</td></tr><tr><td>Modulus of eigenvalues of [P]</td><td>3.83<br>0.95<br>0.95<br>0.95<br>0.95</td><td>4.02<br>0.77<br>0.52<br>0.52<br>0.03</td><td>12.08<br>0.99<br>0.99<br>0.99<br>0.25</td><td>12.15<br>0.69<br>0.43<br>0.01</td><td>0.93<br>0.93<br>0.58<br>0.58</td><td>0.93<br>0.93<br>0.58<br>0.58</td></tr><tr><td>Stability</td><td>Unstable</td><td>Unstable</td><td>Unstable</td><td>Unstable</td><td>Stable</td><td>Stable</td></tr></tbody></table>

† Impacts/cycle is implied.

‡ $X_{max}$ is the maximum deflection of the primary system with impacts.
<figure><img src="image_5.png" /><figcaption>Figure 4. Comparison of two unequispaced impacts/cycle theoretical results with those of Sadek [6]. Both data sets were computed without regard for stability. —, Sadek's results; ---, particular results of general theory. $ξ=0.004$, $μ=0.0484$, $R_1=R_2=0.6$.</figcaption></figure>

problems. Previous problems involve mainly two equispaced† impacts per cycle of the sinusoidal load, $F_0 \sin \Omega t$ [2–5, 7, 8, 10, 11, 14, 22, 23]. Information pertaining to two unequispaced and three impacts per cycle is very sparse [6, 13]. The limitations are reflected in the comparisons where earlier data is restricted to $d = 37(F_0/K)$ in Figure 3 and Table 1 or, regardless of stability, to two unequispaced impacts/cycle in Figure 4. More comprehensive calculations based upon the general theory are presented conventionally in Figure 3 as continuous curves relating values of $dK/F_0$ to the frequency ratio $r$ or $\Omega/\omega$. Particulars of the parameters employed are given at the top of each figure or in Table 1. In addition the periodicity, location, number and distribution of impacts were presupposed from values assumed by previous authors. Consequently, two values of $\tau$ were considered for the equispaced impact cases of Table 1 [13]. Stability zones displayed in Figure 3 are similar to ones measured for the analogous impact-pair. Associated boundaries are indicated by two curves having a common symbol. Any point within compatible boundaries corresponds to a stable motion with the same periodic form. Further comparisons are possible for two equispaced impacts/cycle alone but they do not add to conclusions from present communications [25].

Table 1 suggests that agreement is normally closest with Masri's results for two unequispaced impacts/cycle. Then the otherwise strong correlations intimate that the single major change in $X_{max}$ from 13.79 to 14.77 arises most likely from a printing error in reference [13]. Similar conformity is demonstrated initially for the equispaced impacts but significant differences seem to develop in the eigenvalues. These differences have

<footnote>† Equispaced indicates that the durations between any three consecutive impacts are identical. More generally, unequispaced impacts where consecutive durations are unequal will be implied when the distribution of impacts is not stated. The future abbreviation of cycle will infer cycle of the sinusoidal load.</footnote>
been shown [25] to stem largely from a slight slip in former manipulations of the
simultaneous equations which affects only the equispaced case. The slip also, plausibly,
explains the premature termination of Masri's equispaced results at r equal to 1.04 on
the line d = 37(F₀/K) in Figure 3. It does not account on the other hand for Masri's not
obtaining the very narrow two unequispaced impacts/cycle stability zone near r equal
to 1.00. The maximum eigenvalue modulus is near unity in this zone. Consequently the
"all-or-nothing" stability criterion could easily produce opposite conclusions from even
minor variations in numerical data. Such variations could well arise from Masri's,
apparently, using single rather than, presumably more accurate, double precision arith-
metic in numerical computations [7]. Apart from these discrepancies, agreement is
reasonable overall in Figure 3 and also at point "a" corresponding to the lone three
impacts/cycle example. The absence in this figure of stable equispaced impacts between
0.960 < r < 0.985 supports the previous observation that unevenly spaced impacts are
most likely when r is near one.

An additional comparison is given in Figure 4 of the results from the degeneration of
the general theory into two unequispaced impacts/cycle and those from a Fourier series
approach proposed by Sadek [6]. Sadek disregarded the question of stability so that the
information is presented without such a determination. Figure 4(a) indicates that values
of dK/F₀ correlate quite closely whereas large differences occur generally in Figure 4(b)
for sin τ at frequency ratios, r, near unity. A detailed examination of Sadek's theory
indicates that the discontinuity in the velocity of the primary mass at an impact is treated
improperly. Consequently, equation (13) of reference [6] implies that the phase angle τ
is independent of the amplitude F₀ of the external force. This assertion is refuted clearly
by Sadek's own results presented as the solid curves in Figure 4. An independent increase
in F₀ reduces dK/F₀ but r remains unaltered at a value, for example, of 0.999. Con-
sequently point B is transposed leftward along the solid curve r = 0.999 to position B',
say, if the periodic motion, despite a change in the ratio of durations, is still to have two
unequispaced impacts/cycle. Points B and B' respectively correspond to points A and
A' in Figure 4(b) which obviously relate to different τ.

In summary, the theory's conformity with previously meagre results is promising but
not convincing. More complete verification is needed. Therefore a mechanical system
was built to simulate the mathematical model as far as possible in order to provide a
further independent test.

3. EXPERIMENTAL DETAILS

The experimental analog of the ideal primary system in Figure 1 was composed of
items 4, 6 and 18 in Figure 5. Slotted mass 6 with adapters, 4, was supported on two
sides by flexible, spring steel strips interposed by a very stiff but light hollow beam, 18.
These supports were made so long, 12.00±0.08 in, that their ends closest to the slotted
mass never moved more than 0.38°±0.05° from the vertical. Consequently gravitational
effects should be negligible. The supports were also much stiffer in torsion than flexure
so that the experimental movement was almost uniaxial to conform to conditions imposed
by theory. Tightened screws and bolts were used to connect structural components and
enable the primary system to move integrally. It is practically impossible to prevent all
relative motions, however, so that friction developed at joints. Dissipation was found
from conventional free decay and sinusoidal resonance tests [24, 30] to be equivalent
to a viscous damping ratio, ξ, of 0.0114±0.0005. In addition the first and second natural
frequencies of the primary system were found to be 19.87±0.03 Hz and 430±2 Hz.
Further measurements of select amplitudes and phases of the displacement of the primary
<figure><img src="image_6.png" /><figcaption>Figure 5. Details of experimental apparatus and instrumentation. ① Electromagnetic shaker, Ling Model 400 Series 192; ② Hewlett-Packard function generator 3110B; ③ amplifier Ling Model PA300; ④ aluminium adapter; ⑤ secondary mass made from stainless steel; ⑥ slotted primary mass of die tool steel; ⑦ Brüel & Kjaer accelerometer Type 4333; ⑧ Wayne-Kerr capacitance transducer; ⑨ support for capacitance probe; ⑩ two light cotton threads; ⑪ solid frame; ⑫ Wayne-Kerr feed back amplifier TE MK II; ⑬ Hewlett-Packard multimeter 5306A; ⑭ Tektronix double beam storage oscilloscope 7313; ⑮ adjusting mechanism; ⑯ cold rolled plate; ⑰ heavy concrete block; ⑱ hollow beams; ⑲ weak spring; ⑳ Brüel & Kjaer impedance head, type 8000; ㉑ Brüel & Kjaer conditioning amplifier type 2626; ㉒ Brüel & Kjaer conditioning amplifier type 2626; ㉓ Brüel & Kjaer mass compensation unit type 5565; ㉔ Hewlett-Packard multimeter 5306A.</figcaption></figure>

mass also indicated good linearity over extreme load characteristics. Hence, the primary system fairly represented a linear oscillator for the frequencies between 10 and 30 Hz employed subsequently. The total weight was equivalent to the tip weights of brackets, slotted mass and accelerometer, labelled 7 in Figure 5, plus 33/140 times the weight of each support for an essentially fundamental cantilevered action [24]. This total was calculated from straightforward measurements to be 1.61 lb(f). The corresponding stiffness was approximated for light damping as the square of the measured fundamental radian frequency multiplied by the equivalent weight divided by the gravitational constant. This procedure gave K in Figure 1 as 65.01 lb(f)/in.

The secondary system consisted of a dumbbell-like mass, item 5 in Figure 5, strung by two light threads, 10, from a solid frame. It was identical to the mature stainless steel secondary mass described previously in similar impact-pair experiments [1]. The dumbbell reduced extraneous lifting which happened increasingly at higher loading amplitudes but the threads alleviated rotations about the vertical. Hence, the theoretically assumed free, uniaxial motion was duplicated reasonably for the secondary mass which weighed 0.044 ± 0.002 lb(f). This mass was effectively rigid too because its lowest natural frequency was 500 times, or more, greater than any loading frequency. High speed photography measurements of coefficients of restitution have been reported earlier [1, 25]. Experiments were performed only with magnetically insensitive, mature steel which gave the most consistent results for impact-pairs. As before, unequal coefficients of restitution were achieved by taping one side of the primary mass.

The primary system was connected as in Figure 5 through an impedance head and spring to the electromagnetic shaker. The spring was much weaker than the primary system's effective stiffness so that a force rather than a displacement-like input was obtained. The sinusoidally time varying input, generated to within 2% of the nominal amplitude, was monitored with the aid of the directly coupled impedance head. Compensation for the extraneous mass of this head was achieved by employing the electrical compensation circuit recommended by its manufacturer [31]. Measurements were taken
at night to avoid extraneous building vibrations. Absolute displacements and accelerations of the slotted primary mass were observed separately by using a capacitance displacement transducer, 8, and accelerometer, 7. Displacements were determined relative to one of the aluminium brackets, 4. The other bracket maintained symmetry and prevented unbalance. Experimental procedures have been detailed previously [1, 25] and, hence, will not be described again. The only additional limitation was imposed by the restricted force of the shaker [32] which curtailed measurements at the highest periodic impact numbers. Stability had to be demonstrated conclusively again within an arbitrarily imposed but previously satisfactory time limit of 2 min [1].

4. COMPARISON OF EXPERIMENTAL AND THEORETICAL RESULTS

Experimental and theoretical stability zones are compared in Figures 6–8 with the inverse of $dK/F_0$ in Figure 3 taken as the ordinate. Vertical lines correspond invariably to experimental data whilst the continuous curves present related theoretical results. Coefficients of restitution both average 0.75 in Figures 6 and 7 whereas Figure 8 shows the results of taping one contact surface to produce dissimilar coefficients of 0.75 and 0.61. In actuality, a coefficient of restitution was not constant but decreased generally by about 9% as the frequency ratio, $r$, increased within the experimental range indicated in the abscissa. Variations from the constant coefficients assumed can be anticipated from analogous impact-pairs to give various but on the whole small fluctuations for the different stability zones considered [1, 25]. Other parameters observed experimentally and employed in the numerical computations are reported in the figure captions. Figures 6 and 7 relate to the same vibroimpact system; the two impacts/cycle stability zones were separated only to improve clarity. The experimental five impacts/cycle stability zone, distinguished in Figure 7 by vertical lines terminated with dots, was abbreviated by the restricted force capability of the shaker. Corresponding theoretical extremities are indicated simply by dots. Periodicity, impact rate and temporal distribution needed

<figure><img src="image_7.png" /><figcaption>Figure 6. Comparison of experimental and theoretical stability zones for two impacts/cycle with identical coefficients of restitution. Gap size is 0.0105 ± 0.001 in. (a) Equispaced impacts; (b) unequispaced impacts. |, Range of experimental results. Values used theoretically: $\xi$=0.0114, $\mu$=0.028, $R_1$=$R_2$=0.75. Symbols ○, □ as in Figure 3.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Figure 7. Comparison of experimental and theoretical stability zones for three, four and five impacts/cycle with identical coefficients of restitution. Gap size is 0.0105±0.001 in. |, Range of experimental results. Values used theoretically: $\xi=0.0114$, $\mu=0.028$, $R_1=R_2=0.75$. ●●, Theoretical stability boundary; ●—●, experimental stability zone. Other symbols as in Figure 3.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Figure 8. Comparison of (a) experimental and theoretical stability zones for two and three impacts/cycle and (b) experimental stability zone of sliding-like phenomenon. Coefficients of restitution are unequal and the gap size is 0.0081±0.001 in. |, Range of experimental results. Values used theoretically: $\xi=0.0114$, $\mu=0.028$ and $R_1=0.75$, $R_2=0.61$. Symbols as in Figure 3.</figcaption></figure>
<figure><img src="image_10.png" /></figure>

Figure 9. Typical experimental motions of primary mass with unequal coefficients of restitution. Top and bottom trace in a figure represent the external load and absolute acceleration of the primary mass, respectively. (a) Two unequispaced impacts/cycle with one impact on steel and one on tape; (b) three impacts/cycle with one impact on steel and two on tape; (c) one impact on steel and many on tape; (d) two impacts on steel and many on tape. Horizontal time scale is 16·7 ms/division.

at the start of calculations were obtained from experimental histories like those displayed in Figure 9. (The parts (a)–(c) of Figure 9 correspond to the identically lettered points in Figure 8.) A collision is identified by a sudden acceleration in Figure 9 followed by a gradually decaying envelope as the resulting elastic stress waves diminish. These later transients appeared unimportant providing they were attenuated completely before the subsequent collision [33]. Then the idealization of effectively an acceleration pulse with infinitesimal duration is reasonable. Furthermore, the number and timing of impacts in, for example, Figures 9(a) and (b) can be determined straightforwardly. Difficulties arose, however, when collisions occurred repeatedly and very rapidly on the taped side of the primary mass. This phenomenon is illustrated in Figures 9(c) and (d). It is analogous to the sliding observed for the impact-pair with unequal coefficients of restitution [1]. Such behaviour is outside the scope of the present theory.

Figures 6–8 indicate that the non-sliding experimental and theoretical stability zones generally agree well. The largest differences occur in the highest numbered zones and, somewhat surprisingly, in the two equispaced impacts/cycle zone where coefficients of restitution are equal. A close inspection of the equispaced impact zone revealed that the maximum absolute eigenvalue of $[P]$ in equation (32) varied between 0.996 and 0.999. Computations were repeated with single precision arithmetic for this case alone and the ensuing upper stability boundary can be seen from Figure 6(a) to be clearly different. Correlation with experiment is usually poorer with single precision although double precision arithmetic gives an apparently superfluous spike in the practically significant region of $r$ equal to unity. It was impossible, however, to control and load sufficiently finely there to obtain reliable experimental results which might confirm the spike's existence. A similar difficulty was experienced previously [4] which suggests that this point is worthy of further investigation. These experiences together with the narrowness and closeness of all stability zones in the peaks or troughs in Figures 3 and 6–8 indicate that the system's behaviour in the vicinity of $r$ equal to one is very sensitive to small fluctuations in the parameters. Consequently, the traditional two equispaced
impacts/cycle representation would appear to have limited applicability in such circumstances. Remaining differences between experiment and theory are similar to those noted previously for an impact-pair [1]. Consequently, sources of experimental error introduced by the new spring coupling to the shaker and interactions between closer modes of the primary system seem acceptable. Major discrepancies at the largest loads which excite the highest impact rates are largely created again by the uncontrolled lifting of the experimental secondary mass.

Information presented virtually conventionally in Figures 6–8 cannot be applied easily to what is, probably, the most important practical application where vibroimpacts act as a “damper”. Then, the secondary mass should oppose and, hence, reduce the originally resonant displacement of the primary system. Attenuation happens only when $X_{max}/A$, the ratio of the maximum deflections of the primary system with and without impacts, is less than unity. Consequently the alternative representation of Figure 10, where

<figure><img src="image_11.png" /><figcaption>Figure 10. Comparison of the theoretical and experimental displacement performance of an impact damper with (a) identical coefficients and (b) different coefficients of restitution. Values used theoretically: $\varepsilon=0.0114$, $\mu=0.028$. (a) $R_1=R_2=0.75$; (b) $R_1=0.61$, $R_2=0.75$. Symbols as in Figure 3.</figcaption></figure>

$F_0/(dK)$ is plotted for a constant $r$ near unity against $X_{max}/A$ rather than a variable $r$, is most appropriate. The figure suggests that the primary system's maximum periodic deflection is attenuated most when $r$ equals unity. It then seems beneficial to have equal coefficients of restitution and two stable equispaced impacts/cycle. The attenuation at a given $F_0/(dK)$ diminishes for a small increase in $r$ from unity but, very importantly, amplifications always occur and may be fairly large for a slight decrease in $r$. Therefore the performance of a traditionally designed, virtually viscous-less impact damper [7, 21] will deteriorate noticeably if $r$ falls consistently below one in practice. Consequently the driving frequency or speed will require close controls with the addition of an impact damper.
## 5. CONCLUSIONS

A more general theory has been developed to accommodate any number of repetitive impacts influencing a, formerly, resonant system. Although most incongruences with earlier, sparse predictions of stability may be rationalized, a completely convincing case is not presented for the theory's accuracy. Credibility is enhanced greatly therefore by the close correlation demonstrated between the results from the general theory and comprehensive experiments. The recommended procedure for designing impact dampers is seen to be reasonable for one example of a lightly damped primary system with external sinusoidal loading. This particular example, however, suggests that the operation of the primary system should be maintained at or slightly above its individual fundamental resonance. Careful control of the driving speed may be needed in practice because of the extreme sensitivity displayed by the response in this frequency region.

## ACKNOWLEDGMENTS

The authors gratefully acknowledge the support by the Natural Sciences and Engineering Research Council of Canada. The second author was supported also by an Indian Government Scholarship.

## REFERENCES

1. C. N. BAPAT, N. POPPLEWELL and K. MCLACHLAN 1983 *Journal of Sound and Vibration* **87**, 19-40. Stable periodic motions of an impact-pair.

2. P. LIEBER and D. P. JENSEN 1945 *Transactions of the American Society of Mechanical Engineers* **67**, 523-530. An acceleration damper: development, design and some applications.

3. C. GRUBIN 1956 *Journal of Applied Mechanics, Transactions of the American Society of Mechanical Engineers* **78**, 373-378. On the theory of the acceleration damper.

4. R. N. ARNOLD 1956 *Proceedings of International Congress of Applied Mechanics* **9**, 407-418. Response of an impact vibration absorber to forced vibration.

5. G. B. WARBURTON 1957 *Journal of Applied Mechanics, Transactions of the American Society of Mechanical Engineers* **79**, 322-324. On the theory of the acceleration damper.

6. M. M. SADEK 1965-1966 *Proceedings of the Institution of Mechanical Engineers* **180**, 895-906. The behaviour of the impact damper.

7. S. F. MASRI 1965 Ph.D. Thesis, California Institute of Technology. Analytical and experimental studies of impact dampers.

8. S. F. MASRI and T. K. CAUGHEY 1966 *Journal of Applied Mechanics, Transactions of the American Society of Mechanical Engineers* **33**, 586-592. On the stability of the impact damper.

9. H. DITTRICH 1966 *Ingenieur-Archiv* **35**, 150-171. Study of a discontinuous working vibration damper (in German).

10. D. M. EGLE 1967 *Journal of Engineering for Industry, Transactions of the American Society of Mechanical Engineers* **89**, 653-661. An investigation of an impact vibration absorber.

11. M. I. FEIGIN 1967 *Journal of Applied Mathematics and Mechanics* **31**, 166-171. Dynamic theory of a controlled vibration damping model.

12. S. F. MASRI 1968 *Journal of the Acoustical Society of America* **45**, 1111-1117. Analytical and experimental studies of multi-unit impact dampers.

13. S. F. MASRI 1970 *Journal of the Acoustical Society of America* **47**, 229-237. General motion of impact dampers.

14. V. I. BABITSKII and M. Z. KOLOVSKII 1970 *Mashinovedenie* **2**, 16-24. On the dynamics of a system with an impact vibration absorber (in Russian).

15. S. F. MASRI 1973 *Journal of the Engineering Mechanics Division, Proceedings of the American Society of Civil Engineers* **99**, 669-683. Forced vibration of class of nonlinear dissipative beams.

16. W. M. MANSOUR and D. R. TEIXERIA FILHO 1974 *Journal of Sound and Vibration* **33**, 247-265. Impact dampers with Coulomb friction.
17. C. CEMPEL 1974 *Journal of Sound and Vibration* **34**, 199-209. The multi-unit impact damper: equivalent continuous force approach.

18. S. F. MASRI and K. KAHYAI 1974 *International Journal of Non-linear Mechanics* **9**, 451-462. Steady-state motion of a plate with a discontinuous mass.

19. C. CEMPEL 1975 *Journal of Sound and Vibration* **40**, 249-266. Receptance model of the multi-unit vibration impact neutralizer—"MUVIN".

20. R. K. ROY, R. D. ROCKE and J. E. FOSTER 1975 *Journal of Engineering for Industry, Transactions of the American Society of Mechanical Engineers* **97**, 1317-1324. The application of impact dampers to continuous systems.

21. E. SKIPOR and L. J. BAIN 1980 *Journal of Mechanical Design, Transactions of the American Society of Mechanical Engineers* **102**, 338-343. Application of impact damping to rotary printing equipment.

22. A. YE. KOBRINSKII 1969 NASA Technical Translation TTF-534. Mechanisms with elastic couplings, dynamics and stability.

23. A. E. KOBRINSKII and A. KOBRINSKII 1973 *Vibroimpact System*. Moscow: Nauka Press (in Russian).

24. W. T. THOMSON 1965 *Vibration Theory and Applications*. New Jersey: Prentice-Hall, Inc.

25. C. N. BAPAT 1981 Ph.D. Thesis, University of Manitoba. A study of vibroimpact systems.

26. *IMSL Library* 0007 1979 International Mathematics and Statistical Library, Inc.

27. C. C. FU and B. PAUL 1968 *International Journal of Solids and Structures* **4**, 897-905. Stability of motion of impact tools.

28. C. C. FU 1969 *Journal of Applied Mechanics, Transactions of the American Society of Mechanical Engineers* **91**, 743-749. Dynamic stability of an impact system connected with rock drilling.

29. I. I. BYKHOVSKY 1980 *Fundamentals of Vibration Engineering*. New York: Robert E. Krieger Publishing Company.

30. C. C. KENNEDY and C. D. P. PANCU 1947 *Journal of the Aeronautical Sciences* **14**, 603-625. Use of vectors in vibration measurement and analysis.

31. *Instructions and applications manual for Impedance heads Type 8000 and 8001, 1971 Copenhagen: Brüel & Kjaer.*

32. *Instruction manual for 400 series vibrator and amplifier model PA300. 1972 Hertfordshire, England: Ling Dynamic Systems Limited.*

33. Y. MUZYKA 1981 B.Sc. (*Engineering*) Thesis, University of Manitoba. An experimental investigation of multiple unit impact dampers.

## APPENDIX

Coefficients in the 2N linear simultaneous equations (26) and (30) and elements of the component stability matrix $[P_{ia}]$ in equation (32) are given explicitly here. The following variables are needed in addition to the $C_{1i}$ and $C_{2i}$, $i = 1, 2, ..., N$, defined in equations (9) and (10) of the main text:

$$G_{1i} = \frac{(1 + \mu)R_i\Omega}{(1 + R_i)(\alpha_{(i+N)} - \alpha_{(i+N-1)})},$$

$$G_{2i} = \frac{(1 - \mu R_i)\Omega}{(1 + R_i)(\alpha_{(i+1)} - \alpha_i)},$$

$$G_{3i} = \frac{(R_{(i+1)} - \mu)\Omega}{(1 + R_{(i+1)})(\alpha_{(i+1)} - \alpha_i)},$$

$$G_{4i} = \frac{(1 + \mu)\Omega}{(1 + R_{(i+1)})(\alpha_{(i+2)} - \alpha_{(i+1)})}, \quad i = 1, 2, \dots, N. \quad (A1)$$

Then, for $i = 1, 2, ..., N$,

$$W_{1i} = (C_{1i}/\eta\omega)(G_{2i} - G_{1i} - \xi\omega) - C_{2i},$$

$$W_{2i} = 1 - C_{1i}G_{2i}/\eta\omega,$$

$$W_{3i} = C_{1i}G_{1i}/\eta\omega,$$

$$W_{4i} = C_{1i}\{(r/\eta) \cos \alpha_i + (\xi/\eta) \sin \alpha_i\} + C_{2i} \sin \alpha_i - \sin \alpha_{(i+1)},$$

$$W_{5i} = C_{1i}\{-(r/\eta) \sin \alpha_i + (\xi/\eta) \cos \alpha_i\} + C_{2i} \cos \alpha_i - \cos \alpha_{(i+1)},$$

$$W_{6i} = (C_{1i}/\eta\omega)(G_{1i}\{Y_{(i+N)} - Y_{(i+N-1)}\} + G_{2i}\{Y_{(i+1)} - Y_i\}), \quad (A2)$$
$$V_{1i} = -G_{3i} - C_{4i} - (C_{3i}/\eta\omega)(G_{1i} - G_{2i} + \xi\omega),$$

$$V_{2i} = G_{3i} - G_{4i} - C_{3i}G_{2i}/\eta\omega, \quad V_{3i} = G_{4i}, \quad V_{4i} = C_{3i}G_{1i}/\eta\omega,$$

$$V_{5i} = C_{3i}\{(r/\eta) \cos \alpha_i + (\xi/\eta) \sin \alpha_i\} + C_{4i} \sin \alpha_i - \Omega \cos \alpha_{(i+1)},$$

$$V_{6i} = C_{3i}\{-(r/\eta) \sin \alpha_i + (\xi/\eta) \cos \alpha_i\} + C_{4i} \cos \alpha_i + \Omega \sin \alpha_{(i+1)},$$

$$V_{7i} = -G_{3i}(Y_{(i+1)} - Y_i) - G_{4i}(Y_{(i+2)} - Y_{(i+1)}) + (C_{3i}/\eta\omega) \times (G_{1i}\{Y_{(i+N)} - Y_{(i+N-1)}\} + G_{2i}\{Y_{(i+1)} - Y_i\}), \quad (A3)$$

$$C_{3i} = \omega(\eta C_{2i} - \xi C_{1i}) \quad \text{and} \quad C_{4i} = -\omega(\eta C_{1i} + \xi C_{2i}). \quad (A4)$$

variables $\xi, \omega, r$ and $\eta$ used in the above expressions are defined in equation (3). The $\alpha_{ib}, \mu$ and $Y_i$ are given, respectively, by equations (7), (14) and (18). Periodicity conditions (24) and (25) are employed to bring subscripts within the range 1, 2, ..., N and then all the coefficients are calculable.

Explicit expressions given by Masri [13] for elements of the component matrix $[P_{ia}]$ in equation (32) were modified slightly to accommodate different coefficients of substitution. For example, coefficients $k_r, r = 1, 2, 3, 4$, in reference [13] have to be generalized to $k_{ri}, r = 1, 2, 3, 4$ and $i = 1, 2, ..., N$, where

$$k_{1i} = (1 - \mu R_{(i+1)})/(1 + \mu), \quad k_{2i} = \mu (1 + R_{(i+1)})/(1 + \mu), \\ k_{3i} = (1 + R_{(i+1)})/(1 + \mu), \quad k_{4i} = (\mu - R_{(i+1)})/(1 + \mu). \quad (A5)$$

The following variables are defined to reduce the length of subsequent expressions:

$$S_{1i} = (\xi C_{1i} \cdot \eta C_{2i})/\eta, \quad S_{2i} = C_{1i}/(\omega\eta), \quad S_{3i} = -\omega C_{1i}/\eta, \quad S_{4i} = C_{3i}/(\omega\eta),$$

$$S_{5i} = (\dot{Z}_{ia}/\Omega) - A \cos[\alpha_{(i+1)} + \tau] - S_{6i},$$

$$S_{6i} = \{b_i S_{3i} + (\dot{X}_{ia} - A\Omega \cos[\alpha_i + \tau])C_{3i}/(\omega\eta)\}/\Omega,$$

$$S_{7i} = S_{1i}/S_{5i}, \quad S_{8i} = S_{2i}/S_{5i}, \quad S_{9i} = -1/S_{5i}, \quad S_{10i} = -(\alpha_{(i+1)} - \alpha_i)/(\Omega S_{5i}),$$

$$S_{11i} = (-S_{6i} - AS_{1i} \cos[\alpha_i + \tau] + A\Omega S_{2i} \sin[\alpha_i + \tau] + \dot{Z}_{ia}/\Omega)/S_{5i},$$

$$S_{12i} = A \cos[\alpha_{(i+1)} + \tau] + S_{6i},$$

$$S_{13i} = \{-\omega b_i C_{3i} + (\eta C_{4i} - \xi C_{3i})(\dot{X}_{ia} - A\Omega \cos[\alpha_i + \tau])\}/(\Omega\eta),$$

$$S_{14i} = S_{13i} - A\Omega \sin[\alpha_{(i+1)} + \tau]. \quad (A6)$$

By using equations (A5) and (A6) the elements of $[P_{ia}]$ can be shown to be

$$P_{ia}(1, 1) = S_{1i} + S_{12i}(S_{7i} + S_{9i}), \quad P_{ia}(1, 2) = S_{2i} + S_{8i}S_{12i}, \quad P_{ia}(1, 3) = S_{10i}S_{12i},$$

$$P_{ia}(1, 4) = -AS_{1i} \cos[\alpha_i + \tau] + A\Omega S_{2i} \sin[\alpha_i + \tau] - S_{6i} + S_{11i}S_{12i},$$

$$P_{ia}(2, 1) = k_{1i}(S_{3i} + S_{7i}S_{14i} + S_{9i}S_{14i}), \quad P_{ia}(2, 2) = k_{1i}(S_{4i} + S_{8i}S_{14i}),$$

$$P_{ia}(2, 3) = k_{1i}S_{10i}S_{14i} + k_{2i},$$

$$P_{ia}(2, 4) = k_{1i}(-AS_{3i} \cos[\alpha_i + \tau] + A\Omega S_{4i} \sin[\alpha_i + \tau] - S_{13i}) + k_{1i}S_{11i}(S_{13i} - A\Omega \sin[\alpha_{(i+1)} + \tau]),$$

$$P_{ia}(3, 1) = k_{3i}(S_{6i} + S_{7i}S_{14i} + S_{9i}S_{14i}), \quad P_{ia}(3, 2) = k_{3i}(S_{4i} + S_{8i}S_{14i}),$$

$$P_{ia}(3, 3) = k_{3i}S_{10i}S_{14i} + k_{4i},$$

$$P_{ia}(3, 4) = k_{3i}(-AS_{3i} \cos[\alpha_i + \tau] + A\Omega S_{4i} \sin[\alpha_i + \tau] - S_{13i} + S_{11i}S_{14i}),$$

$$P_{ia}(4, 1) = S_{7i} + S_{9i}, \quad P_{ia}(4, 2) = S_{8i}, \quad P_{ia}(4, 3) = S_{10i}, \quad P_{ia}(4, 4) = S_{11i}. \quad (A7)$$
Variables can be computed in the same manner as the coefficients of equations (A2) and (A3) with the additional terms $\dot{X}_{ia}$, $\dot{Z}_{ia}$, A and $\tau$ determined by the particular periodic solution under consideration. Subscript *i* in equations (A6) and (A7) varies from 1 to *N*.