## 1. Introduction

The vibration isolation system is used to minimize the motions transmitted to the upper part of a system when the base or support of the system moves^{[1]}. The suspension system of car seats is an example of a vibration isolation system^{[2,3]}.

The author’s previous study proposed a vibration isolation system with nonlinear springs and investigated its characteristics^{[4]}. The proposed nonlinear springs are composed of two symmetrical linear springs. These nonlinear springs are characterized by a much simpler structure compared with other nonlinear springs with different structures, such as disc springs^{[5]}, and by the easy adjustability of the spring’s load capacity. The vibration isolation system was constructed using these nonlinear springs and a linear damper. Fig. 1 shows this vibration isolation system. If the nonlinear relationship between the spring force and displacement can be expressed as a polynomial, the response of this nonlinear vibration isolation system can be approximated. The approximation methods for the solutions of the nonlinear equation of motion include the Harmonic Balance Method^{[6,7]} and methods using a higher-order frequency response function^{[8,9]}. In this study, the response component of the system mass at the excitation frequency was measured while the base of the vibration isolation system was moving harmonically and compared with the numerical analysis result.

## 2. Target System

### 2.1 Target System

The nonlinear vibration isolation system in Fig. 1 consists of a mass, two linear springs, and a linear damper. This system has a mass of *m* = 100 kg, a linear spring constant of *k* = 10,000 N/m, a damping coefficient of *c* = 1,000 Ns/m, an original spring length of *l*_{0} = 0.7 m, and an initial spring angle with the horizontal line θ = 0.7752 rad (44.42°).

The equation of motion of the mass when the base of the vibration isolation system receives excitation is as follows:

where *f* is a function of the spring force, which can be found in the reference^{[4]}. The base displacement is given by $x\left(t\right)=X\mathrm{sin}\omega t$
, and the displacement, $y\left(t\right)$
, of the mass is measured from the equilibrium state. $z=y-x$
, Eq. (1) becomes the following:

The response of the mass can be determined by numerically solving this equation for *z* and adding *x* to it. The built-in function ‘ode45’ of MATLAB was used to solve the equation. The response component in multiple frequencies can be found by taking the fft of $y\left(t\right)$
. The response component at the excitation frequency of $\omega =10$
rad/s, and the amplitude of the base displacement, $X=0.05$
m are shown in Fig. 2. This figure shows that the component at 10 rad/s, which is the excitation frequency of the system, and the component at 20 rad/s, which is its harmonic term, are dominant. The component amplitude at the excitation frequency is 0.049328 m.

## 3. Application of the Analysis Method

The Harmonic Balance Method is widely used to approximate the stable or unstable periodic responses of a nonlinear vibration system. In this method, the steady state response of the nonlinear system that receives a harmonic input is represented by the sum of harmonic terms. In a simple case, the response is expressed as follows:

The amplitude of the harmonic term can be determined by substituting the above equation in the nonlinear equation and equating coefficients of the harmonic terms on both sides.

### 3.1 Using One Mode

In Eq. (2), $f\left(z\right)$ represents the nonlinear spring force. After placing the origin of the coordinate axis at a position of static equilibrium, the spring force around this position can be approximated by the third-order polynomial, $f\left(z\right)={k}_{1}z+{k}_{2}{z}^{2}+{k}_{3}{z}^{3}$ . The values of the coefficients can be determined in the range of $-0.13\le z\le 0.109$ using regression analysis as follows:

$\begin{array}{l}{k}_{1}5031.6\text{\hspace{0.17em}N/m,\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}}{k}_{2}=-21126\text{\hspace{0.17em}}{\text{N/m}}^{2},\\ {k}_{3}=-13715\text{\hspace{0.17em}}{\text{N/m}}^{3},\end{array}$

Fig. 3 shows both the original spring force and the regression analysis result; these two curves are almost identical and cannot be distinguished. Using this result, Eq. (2) can be written as follows:

In the above equation, an unknown phase angle, ϕ was used in the excitation force so that one may obtain a fundamental harmonic response containing a single trigonometric term. Assuming $z\left(t\right)=Z\mathrm{sin}\omega t$ , the above equation can be expressed as follows:

If the trigonometric function relation is used, the above equation becomes

To simplify the above equation, $m{\omega}^{2}X$
is represented by *F* . If the coefficients of the $\mathrm{sin}\omega t$
and $\mathrm{cos}\omega t$
terms on both sides are equated, we obtain the following equation:

If we solve the above equations simultaneously, *Z* = 0.04474 m and ϕ = 2.03356 rad. Because $z\left(t\right)=Z\mathrm{sin}\omega t$
is the response to the excitation force, $m{\omega}^{2}X\mathrm{sin}\left(\omega t+\varphi \right)$
, the response to the excitation force, $m{\omega}^{2}X\mathrm{sin}\omega t$
will be $z\left(t\right)=Z\mathrm{sin}\left(\omega t-\varphi \right)$
. If $y\left(t\right)$
is determined by adding $x\left(t\right)=X\mathrm{sin}\omega t$
to $z\left(t\right)$
, the amplitude of the response component at the excitation frequency becomes 0.050044 m, which has an error of 1.45% when compared with the numerically calculated value, 0.049328 m.

In Eq. (6), the constant term, $\frac{1}{2}{k}_{2}{Z}^{2}$ on the left side and the equal sign are not valid. This constant term was generated in the square term of $z\left(t\right)$ . Therefore, to prevent the occurrence of the constant term, the square term was removed from the regression analysis of the spring force as follows:

As a result, we obtained ${k}_{1}=3874.9\text{\hspace{0.17em}N/m,}{k}_{3}=187464\text{\hspace{0.17em}}{\text{N/m}}^{3}$ . The spring force represented by this equation is compared with the original spring force in Fig. 4. It can be seen that the original spring force is not represented properly when regression analysis was performed with the square term removed. When Eq. (9) was used in Eq. (2), which is the equation of motion, and the Harmonic Balance Method was applied as above, the response component amplitude at the excitation frequency became 0.04485 m; thus, the error with the numerically calculated value increased to 9.08%. Therefore, it can be seen that not using the square term in the regression analysis of the spring force to remove the constant term is not a good method.

### 3.2 Using Two Modes

Two modes were used in the Harmonic Balance Method to improve the accuracy of the solution. In other words, as an approximate solution of Eq. (4), the following equation was used:

After substituting $z\left(t\right)$ in Eq. (4) and rearranging the equation, if the coefficients of the $\mathrm{sin}\omega t$ , $\mathrm{cos}\omega t$ , $\mathrm{sin}2\omega t$ and $\mathrm{cos}2\omega t$ terms on both sides are equated, we obtain the following equations:

It is not easy to solve the above equations through a nonlinear equation for the unknown variables ${A}_{1},{A}_{2},{B}_{2},\varphi $ . This equation was solved as follows. If Eqs. (11) and (12) are combined by ${A}_{1}$ and both sides are squared, we obtain the following two equations, where $F=m{\omega}^{2}X$ .

If the above two equations are added, the variable ϕ is removed and ${A}_{1}$ can be expressed as follows:

where

Eq. (13) can be combined by ${A}_{2}$ and expressed for ${A}_{2}$ as follows:

Similarly, Eq. (14) can be expressed for ${B}_{2}$ as follows:

In this way, we obtained three equations, (17), (19), and (20) for ${A}_{1},{A}_{2},{B}_{2}$ respectively. To solve these equations numerically and simultaneously, new ${A}_{1},{A}_{2},{B}_{2}$ were obtained by substituting the initial values of the unknown variables, ${A}_{1},{A}_{2},{B}_{2}$ , in Eqs. (17), (19), and (20), respectively. Then, the unknown values were determined by substituting these values in the equations, and this process was repeated. For the initial value of ${A}_{1}$ , 0.04474 was used, which is the ${A}_{1}$ value when there is only one mode. The initial values of ${A}_{2}$ and ${B}_{2}$ were set to zero. The above calculation process was repeated 100 times. After 10 iterations, the unknown values converged with almost no change. The converged values were ${A}_{1}=0.044739$ , ${A}_{2}=-0.2601105$ , and ${B}_{2}=0.455308$ . Fig. 5 shows the variation of the unknown value ${A}_{1}$ during the repeated calculations, and the unknown values changed very little after 10 iterations.

When the displacement of the mass $y\left(t\right)$ was determined by adding $x\left(t\right)=X\mathrm{sin}\omega t$ to $z\left(t\right)$ , the response component amplitude at the excitation frequency became 0.05008 m. When this was compared with the numerically calculated value, 0.04933 m, the error became 1.52%, which had increased compared with using one mode. This appears to be because, as with the case of using one mode in the approximate equation, the equals sign in the equation of motion was not valid because a constant term that could not be offset appeared in the equation of motion due to the square term used in the regression analysis of the spring force.

## 4. Conclusion

This study approximately analyzed the response of a nonlinear vibration isolation system composed of a nonlinear spring and a linear damper. The displacement of a mass was approximated using the Harmonic Balance Method, while the base of the vibration isolation system moved harmonically. For this purpose, the nonlinear spring force was represented as a third-order polynomial using regression analysis. The response component at the excitation frequency when the response of mass was represented in one mode was approximated and compared with the numerical calculation result.

When the response was approximated using two modes to improve the accuracy of analysis, the accuracy did not increase, unlike the expectation. This appears to be because the second-order term of the spring force in the equation of motion became a constant term, and because this constant term was not offset, the equals sign of the equation of motion became invalid.