Effect of Parameters on the Design of a Suspension System with Four Oblique Springs

(is paper presents the fundamental static and dynamic characteristics of a suspension system consisting of four linear springs arranged in an X-shaped configuration to achieve geometric nonlinearity. (e particular interest is towards the design of a softening spring geometry realizing a quasi-zero stiffness behaviour at large deflections, and the influence of the system parameters is investigated. (e static performance is studied in terms of the force-deflection curve and the dynamic performance in terms of the frequency response curve. (e softening-hardening behaviour of the suspension leads to a frequency response which bends to the lower frequencies reaching a well-defined minimum. It is found that both the static and dynamic behaviours may be described in terms of a single parameter, and a simple closed-form expression is determined which links the damping in the system to the excitation amplitude to achieve the lowest possible resonance frequency.


Introduction
e ever growing demand to enhance the performance of mechanical systems and structures has recently pushed research efforts towards the exploitation of nonlinear effects rather than to their avoidance [1,2]. New solutions and innovative designs have been investigated recently to this purpose, especially in the field of nonlinear dynamics and vibration control [3,4]. Mathematical approaches have been adopted to investigate the dynamic behaviour of nonlinear oscillators, with specific emphasis to prescribed nonlinear functions of stiffness [5,6] and damping [7,8].
Practical applications of the benefits that nonlinearity can introduce into a mechanical system are reported in the field of energy harvesting from vibrations [9], vibration absorbers [10], shock isolators [11], vibration isolators [12], and elastic systems for potential energy increase [13]. In some cases, a nonlinear stiffness element with a quasi-zero stiffness (QZS) characteristic [2,6,8,12,13] has been proposed to cope with the competing requirements of achieving a high-static stiffness to limit the static deflection and a low-dynamic stiffness to improve the dynamic performance. A cubic stiffness characteristic with hardening behaviour has been commonly reported, and a practical mechanical realization consists of a pair of linear springs located perpendicularly to the direction of motion, which incline as the oscillator moves [2]. A global nonlinear behaviour of hardening type assures stability at the equilibrium configuration, where the QZS effect is often required, and this characteristic is obtained by combining elastic elements with positive and negative stiffness in parallel. e exclusive use of a dominant softening stiffness effect, which can be also practically obtained by using magnets arranged in an attractive configuration [14], has the potential disadvantage to eventually lead to a bistable or snap-through behaviour [2,15] for large displacements or even instabilities, which could be undesired or detrimental in some cases.
In this paper, a softening-hardening behaviour is achieved by simply arranging linear springs in an X-shaped mechanical configuration, realizing a quasi-zero stiffness effect for large deflections. Following the preliminary idea presented in [16], a series of research works have recently addressed its exploitation to develop innovative vibration isolators [17][18][19], and it is believed that such characteristic may be of interest to other applications, so that it is further studied in this paper.
Respect to the short communication presented in [16], where the system was excited by a force on the oscillating mass and no effect of gravity was considered, this paper presents a more attractive engineering configuration where the system is base-excited and the effect of the gravity on the static equilibrium is taken into account. is in turn involves the introduction of more design parameters which may be exploited for a better tuning to the desired application. e contribution of the present paper is thus a step forward towards the practical design of the system. A deeper theoretical insight is undertaken which shed new light on the key parameters of the system. A simple and interesting relation between the excitation amplitude and the system damping is uncovered from the study of the frequency response. e analytical and numerical analysis presented in this paper allows a fundamental understanding of the nonlinear characteristics of the system, underpinning prospective pioneering applications which can be potentially attractive to the engineering community for vibration control issues at large [20].

Force-Deflection Characteristic
e model of the system considered in this work is illustrated in Figure 1(a), and it consists of a mass m suspended by two pairs of linear springs arranged in an oblique geometric configuration inside a casing. e pair of springs at the top (bottom) have stiffness k 1 (k 2 ), and the geometry of the spring arrangement is defined by the dimension a and b, as indicated in the figure. e system is excited at its base (the casing) by an imposed displacement y, as illustrated in Figure 1(b), so that the suspended mass can oscillate inside its casing. A viscous damper c is introduced as a generic dissipative term, which is a practical assumption validated by some recent experimental works on similar suspension configurations [17][18][19]. e relative displacement between the mass and the casing is indicated by z, so that when z = 0 the mass is located at the centre of the casing. e system is subject to the acceleration of gravity g, pointing downwards, as indicated in the figure. Due to symmetry, the mass oscillation is assumed to be constrained to the vertical direction, and when the mass moves, the springs rotate changing their length. It is considered that z 0 is the initial displacement achieved upon assembly, as highlighted in Figure 1(a).

General Characteristic.
e relationships between the mass displacement z and the applied static elastic restoring force f s is given by where σ 1 (σ 2 ) is a natural length factor, which accounts for the top (bottom) springs to be assembled in tension or compression respect to their corresponding natural length. is latter is obtained as the geometric distance between each spring attachments at z � z 0 . e spring natural length factor is less (greater) than 1 when the springs are assembled in tension (compression), and it is equal to 1 when the springs are not deformed upon assembly.
Equation (1) may be conveniently rewritten in nondimensional form as z 0 � (z 0 /a), z � (z/a), and p � (mg/ak 2 ). Equation (2) is a highly nonlinear function of the nondimensional displacement z, and its parameters greatly affect the shape of the force-displacement characteristic. As a qualitative indication of the different shapes that the forcedeflection curve may assume, equation (2) is plotted for different values of the system parameters in Figures 2(a)-2(f ), for displacements values in the range − 1.2 < z < 1.2. In all subfigures, the solid lines correspond to the case with k � σ 1 � σ 2 � b � z 0 � p � 1, i.e., the case where the top springs are horizontal in absence of gravity, all springs have equal stiffness, and the casing is as wide as high. Dashed and dotted lines denote variations on one of the parameters, as indicated in each subfigure.
From Figure 2(a), it can be noted that as the stiffness ratio k increases, i.e., the springs at the top become stiffer compared with those at the bottom, this generates a more symmetric behaviour for displacement values around zero. A similar qualitative trend is observed in Figure 2(b) (Figure 2(c)) when the natural length factor σ 1 (σ 2 ) increases (decreases), as these factors have a strengthening effect respect to a change in the corresponding stiffness. As shown in Figure 2(d), the effect of decreasing the aspect ratio b determines an increase of the nonlinear effect as the springs tend to rotate towards the vertical configuration. e effect 2 Shock and Vibration of the initial deflection z 0 is to fundamentally shift the forcedeflection curve mainly along the horizontal axis in a nonlinear manner, as shown in Figure 2(e), while the effect of the nondimensional weight p is to shift the force-deflection curve vertically, as illustrated in Figure 2(f ). Although the analysis performed above, based on a simple observation of Figure 2, is not exhaustive and a more advanced sensitivity study would be required for a deeper investigation, Figure 2 basically illustrates the qualitative change in the shape of the force-deflection curve, as each parameter is varied respect to a given reference configuration.

Symmetric Characteristic.
In fact, the actual objective of this paper is to consider a particular shape of such static force-deflection curve, which motivates the qualitative analysis performed in the previous section. In particular, the interest is towards the case where the system has a symmetric behaviour around the static equilibrium configuration, and the static equilibrium configuration is achieved at z = 0. Such a behaviour is achieved when the force-deflection curve is an odd function of the nondimensional displacement. To impose this condition, equation (2) is expanded in Taylor series around z = 0, and the zero-and second-order coefficients are set to zero yielding Equations (3) and (4) are solved in terms of the stiffness ratio k and the natural length factor σ 2 yielding the following simple relations among the system parameters: Shock and Vibration 3 Shock and Vibration which can be substituted back into equation (2) resulting in the following expression of the static force: Equation (7) only depends on the aspect ratio b, the nondimensional weight p, the natural length factor σ 1 , and the initial deflection z 0 , and it is plotted in Figures 3(a)-3(d) for different combinations of those four system parameters. In particular, the solid lines in all subfigures represent the case where From Figure 3, it is seen that the force-deflection curve is symmetric, with a general softening effect for large deflections. As the natural length factor σ 1 increases, as in Figure 3(a), the force-deflection curve exhibits a quasi-zero stiffness characteristic for large deflections, eventually leading to potential instability (negative stiffness). e effect of the aspect ratio b is seen in Figure 3(b) and is similar to that in Figure 3(a). e effect of z 0 is illustrated in Figure 3(c), where it is seen that the curves are symmetric with respect to z 0 � 1. e effect of increasing p is seen in Figure 3(d) and has an opposite effect respect to an increase in the parameters σ 1 and b in Figures 3(a) and 3(b), respectively.

Symmetric Characteristic with QZS Behaviour.
With the objective to further reduce the system complexity and the number of independent parameters, it is decided to investigate the specific conditions which yields to a QZS behaviour at large deflections (around z ≈ ± 1), followed by a increase in stiffness. Such a situation would guarantee that the mass oscillation will have a global positive stiffness characteristic.
From a mathematical point of view, the condition above is assured when both the first and second derivative of equation (7) are zero, which are the conditions for the appearance of an inflection point with horizontal tangent. Closed-form relationships among the system parameters to attain such a condition were not possible to be achieved; hence, numerical solutions were sought for some combinations of the system parameters. From a first insight, it can be demonstrated that the second derivative of equation (7), respect to the nondimensional displacement, is only dependent on the aspect ratio b, so that a curve is plotted in Figure 4, to show the relation between b and z for is curve is symmetric with respect to the vertical axis, so that only the values for z > 0 are plotted. It may be seen that for b < 1, the inflection point in the force-deflection curve is around z � 1 and dramatically decreases for b > 1.
e system of equations (z 2 /zz 2 )f s � (z/zz)f s � 0 is then numerically solved in terms of z and b, for fixed values of σ 1 and different combinations of z 0 and p, varying in the range 0.1 < z 0 < 0.9 and 0 < p < 2. Results are presented in terms of contour plots in Figure 5, which may be used for design purposes. Figure 5 consists of a table of six panels organized in three rows and two columns, and it should be read as follows: the plots in each of the three rows refer to a different value of σ 1 , while the first and second column show the solution for the inflection point in the force-displacement curve in terms of z and b, respectively, as a function of z 0 and p. It is noted that the inflection point is achieved for displacement values very close to 1 in all cases, while the aspect ratio generally decreases for increasing values of z 0 . Figure 6 shows the corresponding values of the natural length factor σ 2 and the force f s at the inflection point, obtained from Equations (6) and (7), respectively. For σ 1 < 1, it can be seen how σ 2 increases for increasing values of p and decreasing values of z 0 , while f s at the inflection point fundamentally increases for increasing values of z 0 .
e change in f s for changing values of z 0 is opposite for the case when σ 1 > 1than for the case when σ 1 < 1.
As application examples, two different spring configurations are designed using the aid of Figures 5 and 6. In a first case, it is assumed that σ 1 � 1, z 0 � 0.4, and p � 1.5, which corresponds to the circle () in Figures 5 and 6. In particular, from respectively, it is noted that to achieve such configuration it should be set σ 2 � 0.9 approximately, and the expected value of the force f s at the inflection point is about 2.8. e forcedisplacement curve corresponding to this case is plotted from equation (7) in Figure 7(a) as a thick solid curve (other curves will be described further below).
In a second case, it is assumed that it is of interest to suspend a mass with nondimensional weight p � 1 and to set the nondimensional force f s at the inflection point at about 2. We have different design configurations to achieve such requirements. However, if we want the springs at the bottom to be as unloaded as possible upon initial assembly, we can note that the solution indicated by the square () in Figures 5 and 6 would fit the case and corresponds to about e corresponding force-deflection curve is plotted in Figure 7(b) as a thick solid line (other curves will be described further below). Animation S1in the Supplementary Material illustrates the corresponding physical assembly subject to a quasi-static test. e force-deflection curve corresponding to the diamond () and star () in Figures 5 and 6 is shown in Figures 7(c) and 7(d), respectively, as a thick solid curve (other curves will be described further below). Shock and Vibration 5

Dynamic Analysis
With the aim of incorporating the static force-deflection curve with QZS characteristic at large displacements into the dynamic equation of motion of the system and performing an analytical insight, equation (7) is approximated by a polynomial expression. Unfortunately, the classical Taylor series expansion around the static equilibrium configuration would have the limitation that a high order solution would be needed to fit the QZS behaviour at large deflections [16], i.e., far away from the expansion point of the series. To illustrate this, the results of a 7th order Taylor series expansion are illustrated in Figures 7(a)-7(d) as thin dotted lines and compared with the corresponding curve from equation (7), as thick solid lines. It can be observed that a negative stiffness is predicted for values of displacements greater (smaller) than 1 (− 1). To overcome such limitation, the force-deflection curve in equation (7) is approximated by a polynomial expression with the following conditions: (i) the stiffness for z � 0 is set equal to that from Taylor series; (ii) the stiffness at z � 1 (the approximate location of the inflection point) is set to zero; (iii) the derivative of the stiffness at z � 1 is set to zero; and (iv) the force at z � 1 is set equal to that from equation (7). To satisfy these conditions, the following 7th order polynomial expression is obtained: where  Shock and Vibration Equation (8) is plotted in Figures 7(a)-7(d) as a dashed line, and it can be noted that such approximate expression fits the exact force-deflection curve (solid line) better than that from the Taylor series approximation (dotted line). Furthermore, the presence of a negative stiffness behaviour is avoided for nondimensional displacement values around ± 1. e lower panels in each subfigure show the residual errors of the Taylor series and polynomial fitting approximation.
However, from Figures 7(a)-7(d), it can be seen that when |z| becomes greater than 1, the polynomial fitting error rapidly increases. In the following analysis, it is thus assumed that the approximate working condition of the suspension is such that |z| ≤ 1.3, where it is seen that the residual error of the polynomial fitting approximation is below 10%.

Amplitude-Frequency Equation.
e equation of motion of the system depicted in Figure 1is given by where the static spring restoring force f s (z) is given by equation (1), and the overdots denote differentiation respect to time t. e oscillator is excited at the base by a harmonic displacement y � Y cos(ωt + φ) where Y is the displacement amplitude, ω is the angular frequency, and φ is the phase. By using the approximate expression for the spring restoring force given in equation (8), equation (10) may be conveniently written in nondimensional form as where τ � ω 0 tand Ω � (ω/ω 0 ) are the nondimensional time and frequency, respectively, ω 0 � ������ k 2 k 1 /m, ζ � (c/2mω 0 ) is the damping ratio, Y � (Y/a) is the nondimensional displacement amplitude, primes denote differentiation respect to τ, and e expressions given in equation (12) show that the fifth-and seventh-order coefficients are dependent on the cubic one. Both the statics and dynamics of the system exhibiting QZS behaviour for large deflections (at |z| up to approximately 1.3) may be then described in terms of one key parameter only. e values of care calculated based on equation (12) and shown in the contour plots of Figure 8, in the same way as for Figures 5 and 6.
To solve equation (11) in closed form in terms of the amplitude-frequency equation, it is assumed that the system response is predominately harmonic at the excitation frequency, i.e., z � Z cos(Ωτ), and this is substituted back into equation (11), where a first-order harmonic balance approximation is applied to yield with G � 1 + (3/4)cZ 2 + (5/8)δZ 4 + (35/64)εZ 6 .

Backbone Curve and Effect of Damping.
e analytical amplitude-frequency equation reported in equation (13) can be used to investigate the backbone curve and the effect of damping on the system response.
First, the relation between the resonance peak and the corresponding frequency is obtained by setting the discriminant of equation (13)  Shock and Vibration 7 which can be solved in terms of the damping ratio squared to give

Shock and Vibration
By combining equation (15) and equation (13), the expression for the backbone curve, i.e., the unforced and undamped response of the nonlinear oscillator is obtained as  To qualitatively illustrate the shape of the backbone curve for the system considered in this work, equation (16) is plotted in Figure 9(a) as a thick line, for the different values of c.
It can be noted that the frequency of the backbone curve may exhibit a relative minimum and maximum for certain values of c.
ese are due to the transition among the softening and hardening stiffness characteristic in the spring force-deflection curve, and such effect is better investigated below. Of particular interest is the occurrence of a relative minimum in the backbone curve, as this leads to a minimum stable value for the resonance frequency.
An approximate expression for the displacement amplitude at the relative extremes of the resonance frequency in the backbone curve is obtained by differentiating equation (16) with respect to Z, equating to zero, and solving to give where η � ��������� 4 + 6c + c 2 and C (D) indicates the point where a relative minimum (maximum) occurs. e corresponding frequency is determined by substituting equation (17) back into equation (16) to yield (18) e loci of points C (relative minimum) and D (relative maximum) are indicated in Figure 9(a) as a thin solid and dash-dotted line, respectively, while their frequency and amplitude are plotted in Figures 9(b) and 9(c), respectively, as a function of c.
Form Figure 9, it can be noted that for c > 0, the backbone curve presents a relative maximum frequency at a lower displacement amplitude and a relative minimum frequency at a larger displacement amplitude, so that a hardening behaviour manifests at larger amplitudes. In the range − (2/3) < c < 0, the backbone curve presents only a relative minimum frequency, with hardening behaviour at high frequency. In the range − 3 + � 5 √ < c ≤ − (2/3), the backbone curve presents both a relative maximum and a relative minimum frequency, but now the relative maximum appears at a larger amplitude, so that a softening behaviour manifests at larger displacements, potentially leading to instability. For c � − 3 + � 5 √ , the relative maximum and minimum frequencies coincide and an inflection point is formed. For c < − 3 + � 5 √ , no relative maximum or minimum frequency occurs, and an overall softening behaviour manifests, potentially leading to instability for large motion. e minimum frequency and corresponding displacement amplitude, which are indicated by the circle in Figures 9(b) and 9(c), are given, respectively, by 5 − (1/4) ≈ 0.69 and ��������� � 2 + (2/ � 5 √ ) ≈ 1.7. e values of Ω C and Z C are calculated and shown in the contour plots of Figure 10, in the same way as for Figures 5  and 6. It can be noted that for a reasonable range of parameter values, the peak amplitude of the nondimensional displacement does not exceed the approximate value of 1.3, which validates the polynomial fitting approximation of the force-deflection curve presented above.
It is now possible to investigate the effect of damping on the system response. Of particular interest is the value of damping which leads to the minimum resonance frequency, corresponding to point C, discussed above. To this purpose, the ratio (ζ C /Y) is calculated from equations (15) and (17), for different values of excitation amplitude, and is plotted in the contour plots of Figure 11. In Figure 11, dashed lines correspond to Y � 0.01, solid lines correspond to Y � 0.3, and dash-dotted lines correspond to Y � 0.5. As can be noted, despite the large variation in the amplitude of excitation, the ratio (ζ C /Y) does not change much and stays in the range between about 0.28 and 0.38. A closed-form expression can be obtained by expanding equation (15) in Taylor series for small excitation amplitudes and rearranging to give which can be further expressed for small values of c as is is a very simple closed-form expression and shows that the damping to achieve the minimum resonance frequency is proportional to the amplitude of excitation, through a coefficient which is a linear function of the parameter c.
To validate the approximate expression in equation (20), this is plotted in Figure 12 as a dashed line and compared with the results from equation (15) for a relatively high and low value of excitation amplitude. It can be seen that despite the large variation of the excitation amplitude (10% and 50% of the displacement at the QZS point), the approximate expression captures the main trend for values of c greater that approximately − 0.5.  (13), which is quadratic in Ω 2 , and with damping ratio determined according to equation (20). e stability of the steadystate solutions is calculated by applying Floquet theory as detailed in [21] so that stable solutions are represented by solid lines and unstable solutions are represented by dashed lines in the FRCs of Figure 13. For validation purposes, numerical solutions are indicated by circles, and they are obtained by direct integration of the equation of motion of the system, where the exact expression for the static spring restoring force from equation (2) is used. In the case where a multistable solution is expected, different initial conditions were selected to validate the response on both the lower and higher amplitude branches of the frequency response curve.  It can be noted that the approximate FRCs capture the dynamic behaviour of the system reasonably well. Also, it can be seen that the effect of increasing the amplitude of excitation is to widen the bandwidth around the resonance frequency, but this frequency (and its corresponding displacement amplitude) is held constant through the proper selection of the damping according to equation (20). Animation S2 in the Supplementary Material illustrates the physical assembly of the suspension system corresponding to the case in Figure 13(b), when a sinusoidal base excitation at Ω � 1 is provided.

Conclusions
is paper has investigated the static and dynamic characteristics of a nonlinear suspension consisting of four linear springs arranged in an X-shaped configuration to achieve softening characteristics and a QZS behaviour at large deflections. e analytical insight on the behaviour of the system has allowed to highlight the fundamental design strategy to achieve the desired performance. e global force-deflection curve is approximated by a seventh-order polynomial expression, and it is found that the nondimensional force-deflection characteristic with desired QZS behaviour is only dependent on the cubic stiffness coefficient, i.e., the fifth and seventh stiffness coefficients are expressed in terms of the third one. e approximate forcedeflection curve is incorporated into the equation of motion for dynamic analysis. is is performed analytically in terms of the frequency response of the system, and the effect of the parameters is studied. e investigation into the backbone curve of the harmonic response has highlighted the interesting possibility to tune the damping in the system to achieve the lowest possible resonance frequency. In particular, it is found that such damping is proportional to the excitation frequency and linearly related to the cubic

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e author declares that there are no conflicts of interest regarding the publication of this paper.

Supplementary Materials
Supplementary Materials are provided for this paper as follows. S1: animation showing the static response illustrated in Figure 7(b). S2: animation showing the dynamic response illustrated in Figure 13