## Abstract

In this paper, the forced response of a two degrees-of-freedom (DOF) bilinear oscillator with initial gaps involving inelastic collision is discussed. In particular, a focus is placed upon the experimental verification of the generalized bilinear amplitude approximation (BAA) method, which can be used for the accurate estimation of forced responses for bilinear systems with initial gaps. Both experimental and numerical investigations on the system have been carried out. An experimental setup that is capable of representing the dynamics of a 2DOF oscillator has been developed, and forced response tests have been conducted under swept-sine base excitation for different initial gap sizes. The steady-state response of the system under base excitation was computed by both traditional time integration and BAA. It is shown that the results of experiments and numerical predictions are in good agreement especially at resonance. However, slight differences in the responses obtained from both numerical methods are observed. It was found that the time duration where the DOFs are in contact with each other predicted by BAA is longer than that predicted by time integration. Spectral analyses have also been conducted on both experimental and numerical results. It was observed that in a frequency range where intermittent contact between the masses occurs, super-harmonic components of the excitation frequency are present in the spectra. Moreover, as the initial gap size increases, the frequency band where the super-harmonic components are observed decreases.

## 1 Introduction

In many engineering disciplines, it is becoming important to conduct vibration analysis of damaged structures because the vibration responses can be utilized for structural health monitoring and damage detection. However, in many cases, damaged structures show nonlinearities in their dynamic responses, which hinders the application of traditional vibration tools such as modal analysis and frequency response analysis. In particular, structures with breathing cracks [14] and delaminated composites [5] are known to behave as piecewise-linear (PWL) nonlinear dynamical systems because they involve repetitive opening and closing of mechanical boundaries. Namely, the dynamics of PWL systems include sudden changes in the system characteristics, especially the stiffness of the oscillators. Therefore, even if the PWL system is simple, it shows strong nonlinearities in its dynamic responses.

Since the pioneering work by Shaw and Holmes [6], many studies have been conducted to understand the fundamental nature of PWL systems and to explore their practical importance [711]. For instance, Ing et al. studied the dynamics of PWL systems that are close to glazing incidence [12]. Dyskin et al. studied resonances of impact oscillators [13]. They also studied the dynamics of a chain of bilinear oscillators [10] and showed that there are some patterns of resonance in the chain of oscillators. Andreaus and De Angelis conducted numerical studies on the dynamic response of a single degree-of-freedom (DOF) oscillator excited by a base acceleration and constrained by unilateral constraints [14]. Later they conducted experimental studies on the single DOF oscillator subjected to unilateral constraints under harmonic excitation [15,16]. Furthermore, they investigated the influence of geometrical and mechanical characteristics of isolation and mitigation devices on the vibro-impact system experimentally and numerically [17,18]. Moreover, the computation of nonlinear vibration problems such as contact phenomena using time integration requires a small time-step size to produce accurate results. Recently, efforts for reducing such computational costs have been made by several researchers [1,7,1922]. In particular, an approximation method to produce forced response of bilinear systems has been proposed by Jung et al. [23], which is called the bilinear amplitude approximation (BAA) method. This method was then extended to be able to handle bilinear systems with gaps or prestress between the oscillators [24]. The method is able to produce an accurate approximation of the forced response of bilinear systems, which arise from mechanical systems involving repetitive opening and closing of contacting surfaces, or inelastic collision between the DOF. BAA has also been applied to study the dynamics of high dimensional bladed disks with intermittent contact due to cracks [25,26]. However, the method has never been validated with experimental data. Therefore, in this paper, special attention is given to the experimental verification of the BAA method through comparison between the results obtained by BAA, time integration, and most importantly, experimental data. By providing the experimental verification of BAA, it is expected that it be used for predicting dynamic responses of various bilinear systems with high fidelity. Obtaining high-fidelity experimental data for PWL systems is challenging, mostly because the changes in the dynamics of the PWL systems are non-smooth. Among the PWL systems, several experimental studies have been conducted for impact oscillators. For instance, a forced response study has been conducted on a 1DOF impact oscillator consisting of a mass rolling on a parabolic surface and a rigid stop by Todd and Virgin [27], where chaotic behavior and multiple impacts were observed. Virgin et al. also studied an impact oscillator and showed that rich dynamical behavior can be seen experimentally [28]. Bureau et al. conducted an experimental bifurcation analysis of an impact oscillator [29]. Skurativskyi et al. also investigated a 1DOF impact oscillator with a compliant obstacle numerically and experimentally and reported various rich dynamical behavior including grazing bifurcation and chaos [30]. Numerical and experimental studies have also been conducted for impacting pendulums [31,32]. They investigated the bifurcation behaviors of impacting pendulums from both numerical and experimental studies. To date, however, no experimental work has been conducted on PWL systems involving inelastic collision, which arise in systems involving soft contact, or repetitive opening and closing of mechanical interfaces such as cracks, mechanical interfaces, or delamination. This paper is organized as follows. In Sec. 2, the mathematical model for nonlinear time integration is introduced, and the experimental setup that represents the behavior of the bilinear system is presented. In Sec. 3, results of parameter estimation for the experimental setup are presented. In Sec. 4, the results of numerical and experimental forced responses are presented. In particular, the accuracy and applicability of BAA are discussed in detail. Conclusions of the paper are given in Sec. 5.

## 2 Methodology

This section provides a mathematical model of a 2DOF PWL oscillator, as well as an overview of BAA. The experimental setup for the forced response study is then presented.

### 2.1 Mathematical Model.

First, a mathematical model of the 2DOF PWL system of interest is discussed. The system consists of two independent single DOF systems with an initial gap between the DOFs. This system has two states where the system dynamics is linear, which are referred to as open and closed states. The open state means that the DOFs are not in contact with each other, whereas the closed state means that they are in contact with each other. Figure 1(a) shows the initial state of the 2DOF PWL model of interest where each of the two masses is connected to a base by a spring and a damper. The masses are separated by an initial gap g0. The contact between the masses and the base is frictionless. m1 and m2 denote the masses of the DOFs, c1 and c2 are the damping coefficients, and k1 and k2 are the spring constants. Now, consider that the system is excited by a base excitation x0(t). Then the masses start to move, where x1 and x2 are the displacements of the masses measured from their static equilibria, as is shown in Fig. 1(b). The state of the PWL system switches between the open state and closed state depending on the relative displacement between x1 and x2. Now, consider introducing a gap function, g(x1, x2), that represents the gap between the masses with the existence of an initial gap g0, i.e.,
$g(x1,x2)=g0+x2−x1$
(1)
When g(x1, x2) > 0, the 2DOF system is in its open state with no contact between the oscillators as is shown in Fig. 1(b). On the other hand, when g(x1, x2) ⩽ 0, the system is in the closed state where the oscillators are connected to each other with slight allowance of penetration between the masses, as is shown in Fig. 1(c). Since the system is assumed to follow Newton’s third law, there is a contact force between the masses when they are in the closed state. In this study, the contact force is computed by the penalty method. Namely, the contact force is defined as a product between contact stiffness k′ and the gap function. The equation of motion of the 2DOF PWL system is written as
$m1x¨1+c1x˙1+k1x1=c1x˙0+k1x0+fc(x1,x2)m2x¨2+c2x˙2+k2x2=c2x˙0+k2x0−fc(x1,x2)$
(2)
where fc(x1, x2) is the contact force that hinders the masses from penetration, which is written as
$fc(x1,x2)=k′g(x1,x2)$
(3)
The contact stiffness k′ switches depending on the relative displacement between x1 and x2. Namely, k′ is equal to k* when the oscillator has penetration or k′ = 0, i.e.,
$k′={k*,g(x1,x2)⩽00,g(x1,x2)>0$
(4)
Fig. 1
Fig. 1
Close modal

### 2.2 Bilinear Amplitude Approximation.

The generalized BAA method proposed by Tien and D’Souza [24] is a technique to approximate the forced response of PWL systems efficiently and accurately. The BAA method approximates the nonlinear response of PWL nonlinear oscillators by stitching together the linear responses of the system in its open and sliding states. It is noted that the sliding state in BAA is a generalization of the closed state discussed in Sec. 2.1 for systems with multi-DOFs involving contact where the DOFs can move not only in the direction of contact but also in the directions that are perpendicular to the direction of contact. Contact stiffness is used to minimize the penetration between the contacting DOFs when the system is in the sliding state. For a general n DOF mechanical system with piecewise linearity, the equations of motion of the sliding and open states are written as
$Msx¨s+Csx˙s+Ksxs=F(t)+Fg$
(5)
$Mox¨o+Cox˙o+Koxo=F(t)$
(6)
where xs and xo denote the vectors of generalized coordinates for the sliding and the open systems; Ms, Cs, Ks are the mass, damping, and stiffness matrices for the sliding state; Mo, Co, Ko are the mass, damping, and stiffness matrices for the open state; F(t) is the external forcing vector; and Fg is a constant force caused by the contact stiffness when gaps or prestress exist. Denoting the displacement of the ith DOF measured from the static equilibrium position for the open system as xi and the displacement of the ith DOF measured from the static equilibrium position for the sliding system as $x~i$, a coordinate transformation can be defined, i.e., $xi=x~i+δi$ for i = 1, …, n, where δi is the static equilibrium displacement of the masses caused by the contact stiffness k*. After applying this coordinate transformation to the sliding system only, (5) can be written as
$Msx~¨s+Csx~˙s+Ksx~s=F(t)$
(7)
where $xs=x~s+δ$ and $δ=Ks−1Fg$.
By applying a modal coordinate transformation $x~s=Φsq~s$ and xo = Φoqo on (6) and (7), they can be written as
$ΦsTMsΦsq~¨s+ΦsTCsΦsq~˙s+ΦsTKsΦsq~s=ΦsTF(t)$
(8)
$ΦoTMoΦoq¨o+ΦoTCoΦoq˙o+ΦoTKoΦoqo=ΦoTF(t)$
(9)
where Φs and Φo are the modal matrices for the sliding and open states, respectively, and $q~s$ and qo are the corresponding modal coordinates. Equations (8) and (9) can then be decoupled because of the orthonormality of the modes and are expressed as a set of independent equations as
$q~¨s,i+2ζs,iωs,iq~˙s,i+ωs,i2q~s,i=fs,i(i=1,…,ns)$
(10)
$q¨o,i+2ζo,iωo,iq˙o,i+ωo,i2qo,i=fo,i(i=1,…,no)$
(11)
where ns and no are the number of modes for the sliding and open states, respectively; ζs,i and ζo,i are the modal damping ratios; ωs,i and ωo,i are the undamped natural frequencies; and fs,i and fo,i are the modal forces.
Since Eqs. (10) and (11) are decoupled, the analytical solution of each equation can be expressed as a combination of the linear transient response and the linear steady-state response as follows:
$q~s,i=e−ζs,iωs,it[s1,isinωsd,it+s2,icosωsd,it]+(fs,iωs,i2)cos(ωt−θs,i+α)[1−(ω/ωs,i)2]2+2ζs,i(ω/ωs,i)2qo,i=e−ζo,iωo,it[o1,isinωod,it+o2,icosωod,it]$
(12)
$+(fo,iωo,i2)cos(ωt−θo,i+α)[1−(ω/ωo,i)2]2+2ζo,i(ω/ωo,i)2$
(13)
where s1,i, s2,i, o1,i, and o2,i are scalar coefficients of the linear transient responses to be determined; ωsd,i and ωod,i are the damped natural frequencies corresponding to the natural frequencies ωs,i and ωo,i; and $θs,i=tan−1(2ζs,iωs,iω/(ωs,i2−ω2))$ and $θo,i=tan−1(2ζo,iωo,iω/(ωo,i2−ω2))$. The angle α is the phase difference between the excitation and the steady-state responses. The coefficients α, o1,i, o2,i, s1,i, and s2,i in Eqs. (12) and (13) are the unknown coefficients that are determined by enforcing the compatibility conditions on the physical displacement and velocity at the moment when the system switches its state [24]. The compatibility conditions are enforced by minimizing the differences in displacements and velocities by an optimization solver with respect to those coefficients and the time in the sliding state. This produces an approximate solution of the underlying bilinear system.

The equations of motion in (2) are solved using a traditional Runge–Kutta time integration method and the BAA method introduced in this section, and the results are compared with the measurements in Sec. 4.

### 2.3 Experimental Setup.

This section describes the experimental setup that is capable of simulating and measuring the dynamics of PWL nonlinear oscillators. These experiments involve inelastic collisions of oscillators with an initial gap between the masses. Figure 2 shows the experimental setup. The system consists of two moving components with adjustable masses, which are placed on a linear guide system, and their movements are aligned along a single direction. The entire experimental setup is fixed on an electrodynamic shaker table (SSV-750, San-Esu, Japan) so that the base excitation to the system can be applied. The setup is placed such that the direction of the shaker’s excitation coincides with the direction of the linear guide. Springs are attached to the two sliding components and they are connected to a fixture that is fixed to the shaker table. By adjusting the position where the springs are attached to the fixture, the length of the initial gap between the sliding components can be adjusted. In order to simulate an inelastic collision between the oscillators, viscoelastic rubber sheets (GP-35L, Naigai Rubber Industry, Japan) are attached to the surfaces of the moving components that are in contact with each other during oscillation. By adjusting the mass of the components, the natural frequencies can be adjusted. This is necessary because the masses of the components need to be adjusted so that their natural frequencies differ in order to ensure that the components collide with each other during the oscillation. The displacements of the moving components and the table are measured by laser displacement sensors (IL-300, Keyence, Japan).

Fig. 2
Fig. 2
Close modal

## 3 Parameter Estimation

First, experiments have been conducted to estimate the parameters in (2) by studying the forced response of each single DOF system. This process was necessary to account for the masses of the springs themselves and the viscous frictional forces of the linear guide rail and the ones due to the dynamic deformation of the springs.

First, forced response tests were conducted to obtain forced responses of each single DOF system. Each separate single DOF system was excited by a swept-sine base excitation to find the resonant frequency of the systems described in Sec. 2.3. The responses are shown in Fig. 3. It is noted that these experiments were conducted such that no contact between the masses occurs. Assuming a harmonic base excitation of the form x0 = X0ejωt, the resulting response can also be written as a harmonic function, xi = Xiejωt. By substituting these into (2) with fc(x1, x2) = 0, the displacement amplitude of the single DOF system for the base excitation can be analytically obtained as (see, e.g., Ref. [33])
$Xi=1+(2ζiω/ωi)2[1−(ω/ωi)2]2+(2ζiω/ωi)2X0$
(14)
where ω is the excitation frequency of the shaker, ωi is the natural frequency of the ith oscillator, and ζi is the damping ratio. Based on (14), forced responses were computed with the initial sets of parameters shown in Table 1. These responses are plotted in Fig. 3.
Fig. 3
Fig. 3
Close modal
Table 1

The initial and estimated value of the experimental setup

Symbolsm1(kg)c1(N s/m)k1(N/m)m2(kg)c2(N s/m)k2(N/m)
Initial value1.8425.00386.722.3925.00386.72
Estimated value1.8031.14430.002.4619.13417.78
Symbolsm1(kg)c1(N s/m)k1(N/m)m2(kg)c2(N s/m)k2(N/m)
Initial value1.8425.00386.722.3925.00386.72
Estimated value1.8031.14430.002.4619.13417.78
Second, using the measured and computed responses of the single DOF systems, a weighted nonlinear least-squares method was used to estimate the effective masses, equivalent spring constants, and equivalent damping coefficients of the single DOF systems. The parameter estimation was formulated as an optimization problem that minimized the difference between the measured and computed forced responses at resonance:
$minpifi(pi):=sΔXi2(pi)+rΔωi2(pi)fori=1,2$
(15)
where pi = [mi, ci, ki, s, r], mi is the equivalent mass, ci is the damping coefficient, ki is the spring constant, and s and r are variable weighting coefficients. ΔXi(pi) and Δωi(pi) are defined as
$ΔXi(pi)=Xi(pi)−X~iX~i$
(16)
$Δωi(pi)=ωi(pi)−ω~iω~i$
(17)
where $X~i$ is the measured amplitude of displacement at the resonance and $ω~i$ is the measured resonant frequency. The minimization was then performed based on Eq. (15) with the upper and lower bound constraints as shown in Table 2, and an equality constraint s + r = 1. The minimization was conducted using the function “fmincon” in matlab. The estimated values of the parameters based on the minimization are shown in Table 1, and the calculated forced responses with the estimated parameters are shown in Fig. 3. The weighting coefficients converged to s = 0.76 and r = 0.24 for mass 1 and s = 0.49 and r = 0.51 for mass 2. The interpretation of the estimated values is as follows. First, the optimization process is attempted to increase the values of k1 and k2 or decrease the values of m1 and m2 because the resonant frequencies of the oscillators with the initial parameter values are larger than those obtained from the experiments. Therefore, it is reasonable to see that the values of k1 and k2 are both increased. This also corresponds to the additional stiffness due to the frictional force of the linear guide rail. The increase in m2 makes sense because of the additional inertial effect due to the spring motions. The decrease in m1 was unexpected, however, this is due to the fact that the measured damping of the oscillator 1 was higher than that of the oscillator 2, and the optimization solver is needed to increase c1 and decrease m1 because it needed to increase the value of k1.
Table 2

The constraints for parameter estimation

mi(kg)ci(N s/m)ki(N/m)sr
Upper bound3.0050.00430.0011
Lower bound1.800.00380.0000
mi(kg)ci(N s/m)ki(N/m)sr
Upper bound3.0050.00430.0011
Lower bound1.800.00380.0000

Slight differences between the measured and computed responses with the estimated parameters have been observed. These differences are attributed to the fact that the viscous friction between the linear guide also creates a slight nonlinearity from a frictional effect. However, the predicted amplitudes of displacement at the resonances with the optimized parameter values are in good agreement with the measured amplitudes. This agreement is important because the amplitude of oscillations of the system is used for discussing the validity of BAA.

## 4 Forced Response Analysis

With the parameters obtained in the previous section, nonlinear forced response analyses of the 2DOF system have been conducted using time integration and BAA. The results are then compared with the experimental results.

### 4.1 Convergence Study.

First, convergence studies were conducted to investigate the effects of penalty coefficient k* on the response. The 2DOF system represented by (2) is constructed in the matlab Simulink environment and excited by a swept-sine harmonic excitation with an amplitude of 10 mm and a frequency ranging from 1 to 6 Hz. A fourth-order Runge–Kutta method (ode4) was used as the time integration solver with a fixed time-step size of 0.001 s. In Fig. 4, forced responses of the system are shown for various values of k*. One can see that the value of k* greatly affects the results of the forced response. The results do converge when k* reaches 1.0 × 105 N/m. Therefore, in order to reduce the computational cost while maintaining the numerical accuracy, k* = 1.0 × 105 N/m was used to compute the nonlinear forced responses that are shown in Sec. 4.2.

Fig. 4
Fig. 4
Close modal

### 4.2 Nonlinear Forced Response.

Both experimental forced response tests and numerical forced response calculations have been done for the 2DOF PWL system. In particular, the effects of gap size between the masses on the response and the validity of BAA are discussed. Namely, three gap values are chosen, i.e., g0 = 0 mm, 2 mm, and 11 mm.

First, forced response tests have been conducted with swept-sine base excitation for the 2DOF PWL system. The excitation has been applied for 1–6 Hz.

Next, to obtain the steady-state response of the system, forced response calculations have been conducted by numerical time integration and BAA. At each frequency, the measured displacement at the base is used for evaluating x0 in (2). Also, the estimated parameters discussed in Sec. 3 have been used for the numerical simulations. The contact parameter value of k* for time integration is equal to 1.0 × 105 N/m, which was obtained from the convergence study shown in Sec. 4.1.

First, forced response results for g0 = 0 mm are discussed. Figure 5 shows the displacement amplitude of both masses versus the excitation frequency for g0 = 0 mm. As can be seen, results produced by time integration and BAA are in excellent agreement for the entire frequency range. Moreover, they qualitatively agree with the experimental result especially at resonance, with at most 3% error in the maximum amplitude and 5% in the resonant frequency. However, there are some differences in the results of time integration or BAA and that of the experiment. This may suggest that the mathematical model with the estimated parameter values is not accurate enough to quantitatively capture not only the resonance but also the dynamics with frequencies that are higher or lower than the resonant frequency. To take a closer look at the behavior of the masses at the resonance, the time histories of the gap function for a period of oscillation computed by time integration and BAA are shown in Fig. 6. Note that the horizontal axes of Fig. 6 are normalized with respect to the period of excitation. When the gap function becomes equal to or slightly less than zero, it means that the system is in the closed state. As can be seen, both time integration and BAA produced a gap function between 0 and 10 mm. However, some differences between time integration and BAA are observed. Namely, time integration produced a shorter closed state than that produced by BAA.

Fig. 5
Fig. 5
Close modal
Fig. 6
Fig. 6
Close modal

Second, forced response results for g0 = 2 mm are discussed. Figure 7 shows the displacement amplitude of both masses versus the excitation frequency for g0 = 2 mm. Again, the results computed by time integration and BAA are in excellent agreement for the entire frequency range. Note that the response computed by BAA does not contain results with the excitation frequency below 1.25 Hz or above 5.5 Hz because the contact between the masses do not occur in these frequency ranges. Time histories of the gap function at the resonance are shown in Fig. 8. As can be seen, the duration of the closed state predicted by time integration is still shorter than that predicted by BAA. However, the duration of the closed state predicted by BAA has become slightly shorter than that predicted by BAA for the case with g0 = 0 mm.

Fig. 7
Fig. 7
Close modal
Fig. 8
Fig. 8
Close modal

Finally, the results of the forced responses with a larger gap size, i.e., g0 = 11 mm, are discussed. In Fig. 9, there is a primary resonance at 2.1 Hz and a secondary resonance at 2.7 Hz in the results from experiments. This phenomenon is also observed in time integration with slight differences in the maximum amplitude and frequencies, as seen in Fig. 9. On the other hand, the response produced by BAA is smooth at the resonance and does not show the multiple resonances. Figure 10 shows the time histories of the gap function at 2.4 Hz. The contact duration predicted by time integration is still shorter than that predicted by BAA. However, the discrepancy between the results of time integration and BAA, especially in terms of the duration of the closed state, is smaller than that observed for the cases with g0 = 0 and g0 = 2 mm. This is because the contact duration becomes small as the gap size becomes large, and the contact states predicted by both methods become close to each other. However, the phenomenon of multiple resonances is not observed in the results computed by BAA. This may be due to the fact that BAA does not capture some frequency components in the response, as will be discussed in Sec. 4.3.

Fig. 9
Fig. 9
Close modal
Fig. 10
Fig. 10
Close modal

It is noted that as the initial gap size increases, the difference between the time history obtained by time integration and that obtained by BAA decreases, as can be seen in Figs. 6, 8, and 10. This is because the duration of contact period decreases as the initial gap size increases. As can be seen in Fig. 8(a), the time integration method predicts that the gap function increases up to 4mm after the contact and it stays below 5mm until the normalized time is 0.45. Then it increases again and reaches its maximum. On the other hand, based on the results of BAA as shown in Fig. 8(b), the gap function is forced to be in the closed state when the normalized time is between 0 and 0.4. This strict enforcement of the closed state is achieved by an artificial penalty damping factor that is applied when the masses are in contact with each other [24]. This damping effect is not included in the model used for time integration. Therefore, since this enforcement of the closed state generates the difference between the results obtained by time integration and BAA, if this time period decreases due to an increase in the initial gap, BAA will approximate the dynamics of the system better, as seen in Fig. 10.

### 4.3 Spectral Analysis.

In order to further examine the dynamics of the system, particularly when the masses are in the frequency ranges where intermittent contact takes place, spectral analyses have been applied to both experimental and numerical forced responses. For the experimental data, the wavelet transform has been applied to observe the instantaneous spectrum of the experimental data in frequency because the experiment involved a relatively fast frequency sweep. The bump wavelet was applied as the mother wavelet in the wavelet analysis. For numerical responses, the fast Fourier transform (FFT) was applied to compute the frequency spectrum at each excitation frequency. Note that both the wavelet transform and FFT produce comparable results. Figure 11 shows the spectral distribution of the measured and numerically computed responses of both masses with g0 = 2 mm. For numerically computed responses, both time integration and BAA have been considered. In Figs. 11(a) and 11(d), one can observe an order line that represents the frequency component of the excitation frequency for the entire frequency range. Furthermore, one can also see the super-harmonic components, or the frequency components that are integer multiples of the excitation frequency in the frequency range from 1.5 Hz to 4 Hz for both masses. This means that there is intermittent contact between the masses in that frequency range, because intermittent contact of masses is known to produce super-harmonics. Also in Figs. 11(b) and 11(e), one can see that time integration produces the order line representing the frequency component that is identical to the excitation frequency. Furthermore, one can also observe the super-harmonic components that are observed in the experiments. Note that a subharmonic occurs at about 5 Hz and there is also a zero frequency component present in time integration, but these do not occur in the experiments.

Fig. 11
Fig. 11
Close modal

The results computed by BAA are shown in Figs. 11(c) and 11(f). Again, the spectral distributions computed by BAA do not contain results when the excitation frequency is below 1.5 Hz or above 5.5 Hz, because the contact between the masses do not occur in these frequency ranges. Namely, the response is linear in these frequency ranges, and these linear responses contain the frequency component that is identical to the excitation frequency, which is denoted by the dashed line. As can be seen in Figs. 11(c) and 11(f), in addition to the line corresponding to the excitation frequency, a frequency component that is twice the excitation frequency is also observed. Other than these frequency components, no other super-harmonic components are observed, even when the intermittent contact between the masses is strong. This means that BAA does not produce higher harmonics in the response, even though it captures the behavior of intermittent contact between the masses. This aligns with the key idea behind BAA, which assumes that there are only two states in the dynamics of the bilinear system, i.e., closed and open states. It is interesting to note that BAA does not have any subharmonic or zero frequency components, which, unlike time integration, matches the experiments. Furthermore, it is noted that there is a possibility that higher order nonlinearity has been induced at the contact. For instance, for g0 = 2 mm at 2.4 Hz, the gap function produced by BAA is closed from 0 to approximately 0.4 normalized time and open for the rest of the vibration cycle, as can be seen in Fig. 8(b). This means that the contact was assumed to be perfectly inelastic or soft. On the other hand, from 0 to approximately 0.4 normalized time, the gap function produced by time integration is open with a non-zero gap opening of 3–4 mm. This means that the contact between the masses was not perfectly inelastic. This behavior of the masses generated higher order harmonics in the spectra of the response of both masses. Since we observe the higher order harmonics, both in the results of experiments and of time integration as shown in Fig. 11 such phenomenon likely occurred during the experiments.

Figure 12 shows the spectral distributions of the displacement of both masses with a larger initial gap size, g0 = 11 mm. As seen in Figs. 12(a) and 12(d), measured displacements contain not only the frequency component that is identical to the excitation frequency but also the super-harmonic components especially in the vicinity of the resonances. By comparing Fig. 11(a) or 11(d) with Fig. 12(a) or 12(d), one can see that the width of the frequency bands where the super-harmonic components are observed are narrower in Fig. 12(a) or 12(d), than those observed in Fig. 11(a) or 11(d). This trend is shown even more dramatically in the results for time integration, as shown in Figs. 12(b) and 12(e). Namely, the frequency band where super-harmonic components are observed is much smaller than that observed in Figs. 11(b) and 11(e). Since the occurrence of super-harmonic components in the response indicates that intermittent contact between the masses occurs, this result indicates that when the gap size between the masses increases, the frequency band where intermittent contact occurs decreases. This means that as the gap size increases, the masses tend to make less contact with each other as they vibrate. Indeed, the time history of the gap function at resonance shown in Fig. 10(a) shows that the duration of contact is short, which justifies this argument. Figures 12(c) and 12(f) show the results computed by BAA. As can be seen, the intermittent contact between the masses occurs only in the vicinity of the resonance. This is in agreement with the results given by time integration. This time for BAA, almost no super-harmonic component is observed in the response, which means that even though intermittent contact between the masses occurs in the frequency range, the contact is so mild that it almost has no effect on the dynamics of the masses. This indicates that for the prediction of the dynamics of bilinear systems involving higher harmonic components over a very short time period, such as during impact, further extension of the method may be necessary to fully capture the dynamics. However, it is important to note that BAA does match time integration quite well (see Fig. 9), and there is a bigger difference between the experiment and both computational methods.

Fig. 12
Fig. 12
Close modal

## 5 Conclusions

In this paper, the forced response of bilinear systems with initial gaps involving inelastic collision was investigated by numerical simulation and experiments. A focus was placed upon the experimental verification of the generalized BAA, which was developed for the accurate estimation of forced responses of bilinear systems. By conducting spectral analyses on the experimental results and numerical forced responses of the bilinear system, it was shown that super-harmonic components of the excitation frequency can be observed. Moreover, for a small gap size, time integration appeared to match the higher harmonic content of the experiments better than BAA. This is due to the assumptions in the BAA method that lead to some simplifications in the response that remove some of the higher order harmonics. At the same time, BAA was a better match of the experiments than time integration when considering subharmonics. Furthermore, it was shown that as the gap size becomes large, the frequency range that includes super-harmonic components becomes narrow. To further expand the capability of BAA, extension of the method to account for higher-order dynamics would be beneficial.

## Acknowledgment

This paper is based on work partially supported by the National Science Foundation under Grant No. 1902408, program manager Dr. Robert Landers. Any opinions, findings, and conclusions or recommendations expressed in this paper are those of the authors and do not necessarily reflect the views of the National Science Foundation.

## Conflict of Interest

There are no conflicts of interest.

## Data Availability Statement

The datasets generated and supporting the findings of this article are obtained from the corresponding author upon reasonable request.

## References

1.
Saito
,
A.
,
Castanier
,
M.
, and
Pierre
,
C.
,
2009
, “
Estimation and Veering Analysis of Nonlinear Resonant Frequencies of Cracked Plates
,”
J. Sound Vib.
,
326
(
3–5
), pp.
725
739
.
2.
Andreaus
,
U.
, and
Baragatti
,
P.
,
2011
, “
Cracked Beam Identification by Numerically Analysing the Nonlinear Behaviour of the Harmonically Forced Response
,”
J. Sound Vib.
,
330
(
4
), pp.
721
742
.
3.
Casini
,
P.
, and
Vestroni
,
F.
,
2011
, “
Characterization of Bifurcating Non-Linear Normal Modes in Piecewise Linear Mechanical Systems
,”
Int. J. Non-Linear Mech.
,
46
(
1
), pp.
142
150
.
4.
AL-Shudeifat
,
M. A.
, and
Butcher
,
E. A.
,
2013
, “
On the Dynamics of a Beam With Switching Crack and Damaged Boundaries
,”
J. Vib. Control
,
19
(
1
), pp.
30
46
.
5.
Burlayenko
,
V. N.
, and
,
T.
,
2012
, “
Finite Element Nonlinear Dynamic Analysis of Sandwich Plates With Partially Detached Facesheet and Core
,”
Finite Elements Anal. Des.
,
62
, pp.
49
64
.
6.
Shaw
,
S. W.
, and
Holmes
,
P. J.
,
1983
, “
A Periodically Forced Piecewise Linear Oscillator
,”
J. Sound Vib.
,
90
(
1
), pp.
129
155
.
7.
Jiang
,
D.
,
Pierre
,
C.
, and
Shaw
,
S. W.
,
2004
, “
Large-Amplitude Non-Linear Normal Modes of Piecewise Linear Systems
,”
J. Sound Vib.
,
272
(
3–5
), pp.
869
891
.
8.
Makarenkov
,
O.
, and
Lamb
,
J. S.
,
2012
, “
Dynamics and Bifurcations of Nonsmooth Systems: A Survey
,”
Phys. D Nonlinear Phenomena
,
241
(
22
), pp.
1826
1844
.
9.
Uspensky
,
B.
, and
Avramov
,
K.
,
2014
, “
Nonlinear Modes of Piecewise Linear Systems Under the Action of Periodic Excitation
,”
Nonlinear Dyn.
,
76
(
2
), pp.
1151
1156
.
10.
Dyskin
,
A. V.
,
Pasternak
,
E.
, and
Shufrin
,
I.
,
2014
, “
Structure of Resonances and Formation of Stationary Points in Symmetrical Chains of Bilinear Oscillators
,”
J. Sound Vib.
,
333
(
24
), pp.
6590
6606
.
11.
Saito
,
A.
,
2018
, “
Nonlinear Resonances of Chains of Thin Elastic Beams With Intermittent Contact
,”
ASME J. Comput. Nonlinear Dyn.
,
13
(
8
), p.
081005
.
12.
Ing
,
J.
,
Pavlovskaia
,
E.
, and
Wiercigroch
,
M.
,
2006
, “
Dynamics of a Nearly Symmetrical Piecewise Linear Oscillator Close to Grazing Incidence: Modelling and Experimental Verification
,”
Nonlinear Dyn.
,
46
(
3
), pp.
225
238
.
13.
Dyskin
,
A. V.
,
Pasternak
,
E.
, and
Pelinovsky
,
E.
,
2012
, “
Periodic Motions and Resonances of Impact Oscillators
,”
J. Sound Vib.
,
331
(
12
), pp.
2856
2873
.
14.
Andreaus
,
U.
, and
De Angelis
,
M.
,
2016
, “
Nonlinear Dynamic Response of a Base-Excited SDOF Oscillator With Double-Side Unilateral Constraints
,”
Nonlinear Dyn.
,
84
(
3
), pp.
1447
1467
.
15.
Andreaus
,
U.
,
Baragatti
,
P.
,
De Angelis
,
M.
, and
Perno
,
S.
,
2017
, “
Shaking Table Tests and Numerical Investigation of Two-Sided Damping Constraint for End-Stop Impact Protection
,”
Nonlinear Dyn.
,
90
(
4
), pp.
2387
2421
.
16.
Andreaus
,
U.
,
Baragatti
,
P.
,
De Angelis
,
M.
, and
Perno
,
S.
,
2017
, “
A Preliminary Experimental Study About Two-Sided Impacting SDOF Oscillator Under Harmonic Excitation
,”
ASME J. Comput. Nonlinear Dyn.
,
12
(
6
), p.
061010
.
17.
Andreaus
,
U.
, and
De Angelis
,
M.
,
2019
, “
Experimental and Numerical Dynamic Response of a SDOF Vibro-Impact System With Double Gaps and Bumpers Under Harmonic Excitation
,”
Int. J. Dyn. Control
,
7
(
4
), pp.
1278
1292
.
18.
Andreaus
,
U.
, and
De Angelis
,
M.
,
2020
, “
Influence of the Characteristics of Isolation and Mitigation Devices on the Response of Single-Degree-of-Freedom Vibro-Impact Systems With Two-Sided Bumpers and Gaps Via Shaking Table Tests
,”
Struct. Control Health Monitor.
,
27
(
5
), p.
e2517
.
19.
Saito
,
A.
, and
Epureanu
,
B. I.
,
2011
, “
Bilinear Modal Representations for Reduced-Order Modeling of Localized Piecewise-Linear Oscillators
,”
J. Sound Vib.
,
330
(
14
), pp.
3442
3457
.
20.
Zucca
,
S.
, and
Epureanu
,
B. I.
,
2014
, “
Bi-Linear Reduced-Order Models of Structures With Friction Intermittent Contacts
,”
Nonlinear Dyn.
,
77
(
3
), pp.
1055
1067
.
21.
Zucca
,
S.
, and
Epureanu
,
B. I.
,
2018
, “
Reduced Order Models for Nonlinear Dynamic Analysis of Structures With Intermittent Contacts
,”
J. Vib. Control
,
24
(
12
), pp.
2591
2604
.
22.
Tien
,
M.-H.
, and
D’Souza
,
K.
,
2020
, “
Method for Controlling Vibration by Exploiting Piecewise-Linear Nonlinearity in Energy Harvesters
,”
Proc. R. Soc. A Math. Phys. Eng. Sci.
,
476
, p.
20190491
.
23.
Jung
,
C.
,
D’Souza
,
K.
, and
Epureanu
,
B. I.
,
2014
, “
Nonlinear Amplitude Approximation for Bilinear Systems
,”
J. Sound Vib.
,
333
(
13
), pp.
2909
2919
.
24.
Tien
,
M. H.
, and
D’Souza
,
K.
,
2017
, “
A Generalized Bilinear Amplitude and Frequency Approximation for Piecewise-Linear Nonlinear Systems With Gaps or Prestress
,”
Nonlinear Dyn.
,
88
(
4
), pp.
2403
2416
.
25.
Tien
,
M.
,
Hu
,
T.
, and
D’Souza
,
K.
,
2018
, “
Generalized Bilinear Amplitude Approximation and X–Xr for Modeling Cyclically Symmetric Structures With Cracks
,”
ASME J. Vib. Acoust.
,
140
(
4
), p.
041012
.
26.
Tien
,
M.-H.
,
Hu
,
T.
, and
D’Souza
,
K.
,
2019
, “
Statistical Analysis of the Nonlinear Response of Bladed Disks With Mistuning and Cracks
,”
AIAA J.
,
57
(
11
), pp.
4966
4977
.
27.
Todd
,
M. D.
, and
Virgin
,
L. N.
,
1997
, “
An Experimental Impact Oscillator
,”
Chaos Solitons Fract.
,
8
(
4
), pp.
699
714
.
28.
Virgin
,
L. N.
,
George
,
C.
, and
Kini
,
A.
,
2015
, “
Experiments on a Non-Smoothly-Forced Oscillator
,”
Phys. D: Nonlinear Phenomena
,
313
, pp.
1
10
.
29.
Bureau
,
E.
,
Schilder
,
F.
,
Elmegård
,
M.
,
Santos
,
I. F.
,
Thomsen
,
J. J.
, and
Starke
,
J.
,
2014
, “
Experimental Bifurcation Analysis of an Impact Oscillator—Determining Stability
,”
J. Sound Vib.
,
333
(
21
), pp.
5464
5474
.
30.
Skurativskyi
,
S.
,
Kudra
,
G.
,
Wasilewski
,
G.
, and
Awrejcewicz
,
J.
,
2019
, “
Properties of Impact Events in the Model of Forced Impacting Oscillator: Experimental and Numerical Investigations
,”
Int. J. Non-Linear Mech.
,
113
, pp.
55
61
.
31.
Witelski
,
T.
,
Virgin
,
L.
, and
George
,
C.
,
2014
, “
A Driven System of Impacting Pendulums: Experiments and Simulations
,”
J. Sound Vib.
,
333
(
6
), pp.
1734
1753
.
32.
George
,
C.
,
Virgin
,
L. N.
, and
Witelski
,
T.
,
2016
, “
Experimental Study of Regular and Chaotic Transients in a Non-Smooth System
,”
Int. J. Non-Linear Mech.
,
81
, pp.
55
64
.
33.
Meirovitch
,
L.
,
2001
,
Fundamentals of Vibrations
,
McGraw-Hill
,
New York
.