Dynamical analysis of a forced vibrating planar motion of a spring pendulum

Abstract

This work studies the nonlinear movement of a two degrees-of-freedom (DOF) spring pendulum that is dampened and affected by a harmonic force externally. It is presumed that the spring’s pivot point travels along an elliptic route. Lagrange’s equations are utilized to generate the regulating system of motion. The multiple-scales approach (MSA) is used to gain the system’s analytic solutions up to the third-order approximation. Therefore, all resonance cases that have emerged are categorized, wherein two of them are scrutinized at once. As a result of the removal of secular terms, the solvability constraints are attained and then the steady-state solutions are investigated. The examined motion’s temporal evolution, the resonance response curves, and the solutions at the steady-state are all depicted graphically. In compliance with the Routh–Hurwitz criteria (RHC), all possible fixed points (FPs) for the steady and unsteady cases are found and displayed. The stability zones are examined and analyzed to estimate the effect of various factors on the system’s behavior. This model has gained prominence recently due to its industrial uses in seismic isolation systems for buildings and structures. However, in seismic engineering, a 2DOF vibrating pendulum system can be used as part of a seismic isolation system designed to protect buildings and infrastructure from earthquake-induced vibrations. The pendulum mechanism helps to absorb and dissipate seismic energy, reducing the amount of force transmitted to the structure. During an earthquake, the ground motion acts as an external harmonic force on the building. The distribution of mass and the structural layout can cause rotational moments that act on the building. The pendulum system can be tuned to counteract these moments, helping to stabilize the structure. The pendulum system allows for both horizontal and vertical displacement, providing two degrees of freedom. This capability is essential for accommodating the complex, multi-directional nature of seismic waves.

Introduction

Recently, many studies have appeared that deal with the nonlinear motion of vibratory systems in applied mechanics and study the effects of vibrations on structures, machines, vehicles, etc. The spring pendulum’s dampened motion,^1–6 with linear or nonlinear stiffness, has contributed to addressing some of these problems. The case of internal resonance for the motion of a subjected elastic pendulum to an excitation force along its arm is studied in Reference 1, while the chaotic responses for its motion are examined in References 2 and 3. The positive impact of higher-order extensions on these responses is explored in Reference 4. It is found that the obtained solutions up to the second order outscored on the first approximations when looking at maximum Lyapunov exponents and bifurcation plots. In Reference 5, the procedure of the phase plane and the frequency response equation are used to examine the stability of this system. A modification of this study is examined in Reference 6 when the pendulum’s pivot point has a vertical oscillation.

The approximate solution of a 2DOF auto-parametric system is obtained in References 7 and 8 using the harmonic balancing method (HBM), ⁹ where the chaotic reactions are numerically investigated for this system. Additionally, the investigation looks into the stable and unstable zones. The same methodology utilized was utilized in Reference 10 to investigate how a damped spring system responds to vibration under the impact of a harmonic excitation force. Examining the stability assessments of the obtained results confirms that the numerical outcomes are in high consistency with the approximate one. In Reference 11, the HBM is used to look at the solutions of the equations of motion (EOM) for a severely nonlinear oscillator coupled to a spring’s pendulum. Furthermore, the chaotic behavior of the spring is investigated as well as reducing the vibration near the parametric resonance. The method of averaging ⁹ is utilized in Reference 12 to gain the approximate solution of an oscillating auto-parametric absorber. The limitations of stability and the analysis of bifurcation for the averaging system of EOM are investigated. In Reference 13, the influence of a vibrating auto-parametric absorber is examined. Reference 14 investigates the simple pendulum’s motion connected to a tuned absorber that is excited harmonically. The asymptotic solutions are acquired using MSA ⁹ and represented graphically to show the impact of used parameters on the motion.

It should be highlighted that many academics have become interested in the damped motion of linear or nonlinear elastic pendulums in various pathways.^15–25 In this regard, the planar motion of a harmonic oscillator close to resonance when the pivot point is assumed to move on trajectories of an ellipse and on a closed Lissajous curve is examined in References 26 and 15, respectively. A few limited circumstances for specified motions are provided. In Reference 16, the authors studied the stability analyses of a 3DOF dynamical system of a solid body attached to a damped linear spring with a fixed suspension point. The extensions of this problem are explored in References 17 and 18 with respect to the stiffness of linear and nonlinear springs. The pivot point of the pendulum was constrained to be on an ellipse. In References 19 and 20, a hard body pendulum is further described. Utilizing the MSA, ⁹ the approximations are obtained, and the necessary modulation’s equations are achieved considering the constraints' solvability. The motion of various vibrating systems with the impact of harmonic external torques and forces is explored References 21–25 as well as the nonlinear analysis' stability for these systems.

The dynamical motion of a spring pendulum with weak nonlinearity is discussed in Reference 27, where its suspension point travels in a circular trajectory. The impact of a damper on the spring’s motion is examined in Reference 28 as a general case. The numerical study of this system for periodic, quasi-periodic, and chaotic motion is presented in Reference 29. Whereas in Reference 30 the authors examined the movements of a mass hanging with one end of a nonlinearly damped elastic pendulum under the effect of the spring’s transverse and radial forces, while the other end is considered to be fixed. The improvement of this work was found in Reference 31, where two methods with three time scales are used to decrease the error to a minimum. The EOM’s approximate solutions are obtained using MSA. To fulfill the solvability criteria, the secular term’s elimination condition is preserved. Most of the categories' resonances are studied. The stability criteria of the exact solution of a simple pendulum are explored in Reference 32, in which its fulcrum oscillates both horizontally and vertically. For a few distinct values of the pendulum’s parameters, a comparison of the obtained numerical and analytical solutions (AS) is presented. The periodic motions of a periodically forced nonlinear elastic pendulum are explored in Reference 33 using the semi-analytical approach, as well as the motion’s bifurcation and stability analyses. In Reference ³⁴ , RHC are used to explore the stability and instability zones, which are subsequently studied in terms of steady-state solutions. In Reference 35, the authors obtained analytical solutions of an oscillating 2DOF dynamical system using MSA. All possible FPs are determined in view of the examined resonance cases, while the stability of the triple pendulum is investigated using the nonlinear stability analysis approach in Reference 36. The integrated approach investigates the fractal modification of a nonlinear oscillator relevant to nanotechnology, as explored in prior studies.^37–39 An improved harmonic balance method is used to solve the Duffing-harmonic equation. ⁴⁰ The variational approach is applied to solve different models.^41–44 Some nonlinear oscillators are investigating by using the homotopy perturbation method.^45,46

This article examines the nonlinear behavior of a 2DOF damped spring pendulum when external harmonic forces are applied. The spring’s pivot point is assumed to follow an elliptic path. Lagrangian equations are used to construct the EOM in accordance with the generalized coordinates of the system. Up to the third-order approximation, the desired solutions to these equations are acquired by using the MSA. According to the used perturbation approach, the used parameters may be written in terms of a small parameter. As a result, all discovered cases of resonance are grouped, with two of them examined concurrently. Once the elimination of secular terms, the solvability restrictions are fulfilled, and the solutions at the scenario of steady-state are then analyzed. The steady-state solutions, the resonance response curves, and the time evolution of the investigated motion are all pictorially displayed. All FPs are identified at this case and displayed in accordance with the RHC. To assess the good effects of different parameters on the motion’s behavior, the stability zones are examined and analyzed. Due to its industrial applications, the examined model is now widely used in seismic isolation systems to safeguard buildings and structures during earthquakes. The 2DOF vibrating pendulum system has become a key component in seismic isolation systems, protecting buildings and infrastructure from earthquake-induced vibrations.

Dynamical modeling

In this section, a vibrating pendulum system of a mass $m$ coupled to a damped massless spring with a natural length ${\mathcal{\ell }}_0$ and stiffness $k_1$ and $k_2$ , which are the linear and nonlinear spring stiffness, respectively. Stiffness $k_1$ is typically measured in Newton’s per meter $\left(N.m^{-1}\right)$ , but stiffness $k_2$ is typically measured in Newton’s per cubic meter $\left(N.m^{-3}\right)$ (see Figure 1). The end $N$ is restricted to moving with a regular angular velocity $\Omega$ on an elliptic path with radii $b$ and $a$ . The projections of $N$ on the axes $OY$ and $OX$ are $b\hspace{0.17em}\mathit{sin}\left(\Omega t\right)$ and $a\hspace{0.17em}\mathit{cos}\left(\Omega t\right)$ , respectively. On the dotted circle tangent to the ellipse of radius $b$ , it is noted that the point $N^{\prime }$ corresponds to the projection of the point $N$ on the dotted circle. The planar motion of the system is taken into account as the result of harmonic forces $F\left(t\right)=F\hspace{0.17em}\mathit{cos}\left({\Omega }_{1}t\right)$ , in which ${\Omega }_1$ and $F$ are, respectively, the frequencies and amplitudes of $F\left(t\right)$ . Considering that $x\left(t\right)$ is the dynamic spring’s extension, $c_1$ and $c_2$ are the coefficients of damping in both linear and seesaw vibrations.

Figure 1.

The dynamical model.

With the help of the previous description, the energies of potential $U$ and kinetic $T$ are expressed, respectively, in the form below

$\mathit{U}=\frac{\mathbf{1}}{\mathbf{2}}{\mathit{k}}_{\mathbf{1}}{\left(\mathit{x}+{\mathit{x}}_{\mathit{r}}\right)}^{\mathbf{2}}+\frac{\mathbf{1}}{\mathbf{4}}{\mathit{k}}_{\mathbf{2}}{\left(\mathit{x}+{\mathit{x}}_{\mathit{r}}\right)}^{\mathbf{4}}-\mathit{m}\mathit{g}\left[\mathit{a}\hspace{0.17em}\mathit{cos}\left(\mathbf{\Omega }\mathit{t}\right)+\left(\mathbf{\ell }+\mathit{x}\right)\mathit{cos}\hspace{0.17em}\mathit{\phi }\right],\mathit{T}=\frac{\mathbf{1}}{\mathbf{2}}\mathit{m}{\mathbf{\Omega }}^{\mathbf{2}}[{\mathit{a}}^{\mathbf{2}}\hspace{0.17em}{\mathit{sin}}^{\mathbf{2}}\left(\mathbf{\Omega }\mathit{t}\right)+{\mathit{b}}^{\mathbf{2}}\hspace{0.17em}{\mathit{cos}}^{\mathbf{2}}\left(\mathbf{\Omega }\mathit{t}\right)]+\frac{\mathbf{1}}{\mathbf{2}}\mathit{m}\left[{\dot{\mathit{x}}}^{\mathbf{2}}+{\left(\mathbf{\ell }+\mathit{x}\right)}^{\mathbf{2}}{\dot{\mathit{\phi }}}^{\mathbf{2}}\right]+\mathit{m}\mathbf{\Omega }\dot{\mathit{x}}\left[\mathit{b}\hspace{0.17em}\mathit{cos}\left(\mathbf{\Omega }\mathit{t}\right)\mathit{sin}\hspace{0.17em}\mathit{\phi }-\mathit{a}\hspace{0.17em}\mathit{sin}\left(\mathbf{\Omega }\mathit{t}\right)\mathit{cos}\hspace{0.17em}\mathit{\phi }\right]+\mathit{m}\mathbf{\Omega }\left(\mathbf{\ell }+\mathit{x}\right)\dot{\mathit{\phi }}\times \left[\mathit{b}\hspace{0.17em}\mathit{cos}\left(\mathbf{\Omega }\mathit{t}\right)\mathit{cos}\hspace{0.17em}\mathit{\phi }+\mathit{a}\hspace{0.17em}\mathit{sin}\left(\mathbf{\Omega }\mathit{t}\right)\mathit{sin}\hspace{0.17em}\mathit{\phi }\right].$

(1)

here $\mathcal{\ell }={\mathcal{\ell }}_0+x_r$ in which $x_r$ is the static spring’s protraction, $\phi$ is the rotatory angle, and $g$ indicates the gravity’s acceleration.

The Lagrangian $S=T-U$ of the given model may be easily achieved. Consequently, the following Lagrangian equations can be used to obtain the EOM ²⁹

$\frac{d}{dt}\left(\frac{\partial S}{\partial \dot{x}}\right)-\frac{\partial S}{\partial x}=Q_x,\frac{d}{dt}\left(\frac{\partial S}{\partial \dot{\phi }}\right)-\frac{\partial S}{\partial \phi }=Q_{\phi },$

(2)

where $\left(\dot{x},\dot{\phi }\right)$ and $\left(x,\phi \right)$ represent, respectively, the system’s generalized velocities and coordinates, while $\left(Q_x,Q_{\phi }\right)$ are the generalized forces that can be represented in the following forms

$Q_x=F\hspace{0.17em}\mathit{cos}\left({\Omega }_{1}t\right)\mathit{sin}\hspace{0.17em}\phi -c_1\dot{x},Q_{\phi }=\left(\mathcal{\ell }+x\right)F\hspace{0.17em}\mathit{cos}\left({\Omega }_{1}t\right)\mathit{cos}\hspace{0.17em}\phi -c_2\dot{\phi }.$

(3)

Inserting (1) and (3) in (2), to get the following EOM in a dimensionless form

$\ddot{y}+{\alpha }_1\dot{y}+y+\beta y^3+3\beta y_r{y}^2+3\beta {y_r}^{2}y+{\omega }^2\left(1-\hspace{0.17em}\mathit{cos}\hspace{0.17em}\phi \right)-\left(1+y\right){\dot{\phi }}^2-\frac{{\varpi }^2}{\mathcal{\ell }}\left[b\hspace{0.17em}\mathit{sin}\left(\varpi \tau \right)\mathit{sin}\hspace{0.17em}\phi +a\hspace{0.17em}\mathit{cos}\left(\varpi \tau \right)\mathit{cos}\hspace{0.17em}\phi \right]=f_1\hspace{0.17em}\mathit{cos}\left(p_1\tau \right)\mathit{sin}\hspace{0.17em}\phi ,\ddot{\phi }{\left(1+y\right)}^2+{\alpha }_2\dot{\phi }+2\dot{y}\dot{\phi }\left(1+y\right)+{\omega }^2\left(1+y\right)\mathit{sin}\hspace{0.17em}\phi -\left(1+y\right)\frac{{\varpi }^2}{\mathcal{\ell }}\times [b\hspace{0.17em}\mathit{sin}\left(\varpi \tau \right)\mathit{cos}\hspace{0.17em}\phi -a\hspace{0.17em}\mathit{cos}\left(\varpi \tau \right)\mathit{sin}\hspace{0.17em}\phi ]=\left(1+y\right)f_1\hspace{0.17em}\mathit{cos}\left(p_1\tau \right)\mathit{cos}\hspace{0.17em}\phi ,$

(4)

where the dot indicates the regarding differentiation with $\tau$ and the following dimensionless parameters are furthermore present

$\tau ={\omega }_{1}t,y=\frac{x}{\mathcal{\ell }},y_r=\frac{x_r}{\mathcal{\ell }},{\omega }_1=\sqrt{\frac{k_1}{m}},{\omega }_2=\sqrt{\frac{g}{\mathcal{\ell }}},\omega =\frac{{\omega }_2}{{\omega }_1},\varpi =\frac{\Omega }{{\omega }_1},\beta =\frac{k_2{\mathcal{\ell }}^2}{m{{\omega }_1}^2},{\alpha }_1=\frac{c_1}{m{\omega }_1},{\alpha }_2=\frac{c_2}{m{\mathcal{\ell }}^2{\omega }_1},f_1=\frac{F}{m\mathcal{\ell }{{\omega }_1}^2},p_1=\frac{{\Omega }_1}{{\omega }_1}.$

(5)

The proposed method

Consider checking the oscillations of the system near the position of static equilibrium. Then, functions $\mathit{sin}\hspace{0.17em}\phi$ and $\mathit{cos}\hspace{0.17em}\phi$ may be roughly calculated using Taylor’s series up to third order. Consequently, equations (4) are rewritable in the form

$\ddot{y}+{\alpha }_1\dot{y}+y+\beta y^3+3\beta y_r{y}^2+3\beta {y_r}^{2}y+\frac{{\omega }^2}{2}{\phi }^2-\left(1+y\right){\dot{\phi }}^2-\frac{{\varpi }^2}{\mathcal{\ell }}\left[b\hspace{0.17em}\mathit{sin}\left(\varpi \tau \right)(\phi -\frac{{\phi }^3}{6})+a\hspace{0.17em}\mathit{cos}\left(\varpi \tau \right)(1-\frac{{\phi }^2}{2})\right]=f_1\hspace{0.17em}\mathit{cos}\left(p_1\tau \right)\left(\phi -\frac{{\phi }^3}{6}\right),\ddot{\phi }{\left(1+y\right)}^2+2\dot{y}\dot{\phi }\left(1+y\right)+{\alpha }_2\dot{\phi }+{\omega }^2(1+y)(\phi -\frac{{\phi }^3}{6})-\frac{{\varpi }^2}{\mathcal{\ell }}\left(1+y\right)\times \left[b\hspace{0.17em}\mathit{sin}\left(\varpi \tau \right)(1-\frac{{\phi }^2}{2})-a\hspace{0.17em}\mathit{cos}\left(\varpi \tau \right)(\phi -\frac{{\phi }^3}{6})\right]=\left(1+y\right)f_1\hspace{0.17em}\mathit{cos}\left(p_1\tau \right)\left(1-\frac{{\phi }^2}{2}\right),$

(6)

next, a small parameter $0<\epsilon \ll 1$ is used to represent the coefficients of damping, amplitudes of forces, and the ellipse’s semi axes, as follows

${\alpha }_j={\epsilon }^2{\stackrel{\sim }{\alpha }}_j;j=\left(1,2\right),\beta ={\epsilon }^2\stackrel{\sim }{\beta },f_1={\epsilon }^3{\stackrel{\sim }{f}}_1,a={\epsilon }^2\stackrel{\sim }{a},b={\epsilon }^2\stackrel{\sim }{b}.$

(7)

It is assumed that the oscillations' amplitudes are of order of $\epsilon$ , then we can express the variables $y$ and $\phi$ in the following form

$y=\epsilon z,\phi =\epsilon \varphi .$

(8)

Using MSA, the functions $z$ and $\varphi$ can be addressed as follows ⁹

$z={\displaystyle \sum }_{k=0}^2{\epsilon }^k{z}_{k+1}\left({\tau }_0,{\tau }_1,{\tau }_2\right)+O\left({\epsilon }^3\right),\varphi ={\displaystyle \sum }_{k=0}^2{\epsilon }^k{\varphi }_{k+1}\left({\tau }_0,{\tau }_1,{\tau }_2\right)+O\left({\epsilon }^3\right),$

(9)

where ${\tau }_n={\epsilon }^n\tau \left(n=0,1,2\right)$ performs specific time scales, where ${\tau }_0$ is fast time scales and ${\tau }_1,{\tau }_2$ are slow ones. To formulate $z\left({\tau }_n\right)$ and $\varphi \left({\tau }_n\right)$ , we can write $\tau$ in terms of ${\tau }_n$ as follows ⁴⁷

$\frac{d}{d\tau }=\frac{\partial }{\partial {\tau }_0}+\epsilon \frac{\partial }{\partial {\tau }_1}+{\epsilon }^2\frac{\partial }{\partial {\tau }_2},\frac{d^2}{d{\tau }^2}=\frac{{\partial }^2}{\partial {\tau }_0^2}+2\epsilon \frac{{\partial }^2}{\partial {\tau }_0\partial {\tau }_1}+{\epsilon }^2\left(\frac{{\partial }^2}{\partial {{\tau }_1}^2}+2\frac{{\partial }^2}{\partial {\tau }_0\partial {\tau }_2}\right).$

(10)

It is important to observe that terms of higher order than ${\epsilon }^3$ are omitted for their small size.

Substituting (7)–(10) into EOM (6) and equating the various powers' coefficients of $\epsilon$ with zero to have the following categories of differential equations

Order of ( $\epsilon$ )

$\frac{{\partial }^2{z}_1}{\partial {{\tau }_0}^2}+z_1=0,\frac{{\partial }^2{\varphi }_1}{\partial {{\tau }_0}^2}+{\omega }^2{\varphi }_1=0,$

(11)

Order of $\left({\mathit{\epsilon }}^{\mathbf{2}}\right)$

$\frac{{\partial }^2{z}_2}{\partial {{\tau }_0}^2}+z_2=\frac{\stackrel{\sim }{a}{\varpi }^2}{\mathcal{\ell }}\mathit{cos}\left(\varpi {\tau }_0\right)-\frac{{\omega }^2}{2}{\left({\varphi }_1\right)}^2+{\left(\frac{\partial {\varphi }_1}{\partial {\tau }_0}\right)}^2-2\frac{{\partial }^2{z}_1}{\partial {\tau }_0\partial {\tau }_1},\frac{{\partial }^2{\varphi }_2}{\partial {{\tau }_0}^2}+{\omega }^2{\varphi }_2=\frac{\stackrel{\sim }{b}{\varpi }^2}{\mathcal{\ell }}\mathit{sin}\left(\varpi {\tau }_0\right)-{\omega }^2{z}_1{\varphi }_1-2\left(\frac{\partial z_1}{\partial {\tau }_0}\frac{\partial {\varphi }_1}{\partial {\tau }_0}+\frac{{\partial }^2{\varphi }_1}{\partial {\tau }_0\partial {\tau }_1}+z_1\frac{{\partial }^2{\varphi }_1}{\partial {{\tau }_0}^2}\right),$

(12)

Order of $\left({\mathit{\epsilon }}^{\mathbf{3}}\right)$

$\frac{{\partial }^2{z}_3}{\partial {{\tau }_0}^2}+z_3=z_1{\left(\frac{\partial {\varphi }_1}{\partial {\tau }_0}\right)}^2-2\left(\frac{{\partial }^2{z}_1}{\partial {\tau }_0\partial {\tau }_2}+\frac{{\partial }^2{z}_2}{\partial {\tau }_0\partial {\tau }_1}\right)+2\frac{\partial {\varphi }_1}{\partial {\tau }_0}\left(\frac{\partial {\varphi }_2}{\partial {\tau }_0}+\frac{\partial {\varphi }_1}{\partial {\tau }_1}\right)-{\stackrel{\sim }{\alpha }}_1\frac{\partial z_1}{\partial {\tau }_0}-3\stackrel{\sim }{\beta }{y_r}^2{z}_1-{\omega }^2{\varphi }_1{\varphi }_2+\frac{\stackrel{\sim }{b}{\varpi }^2}{\mathcal{\ell }}\mathit{sin}\left(\varpi {\tau }_0\right){\varphi }_1-\frac{{\partial }^2{z}_1}{\partial {{\tau }_1}^2},\frac{{\partial }^2{\varphi }_3}{\partial {{\tau }_0}^2}+{\omega }^2{\varphi }_3={\stackrel{\sim }{f}}_1\hspace{0.17em}\mathit{cos}\left(p_1{\tau }_0\right)+\frac{\stackrel{\sim }{b}{\varpi }^2}{\mathcal{\ell }}\mathit{sin}\left(\varpi {\tau }_0\right)z_1-\frac{\stackrel{\sim }{a}{\varpi }^2}{\mathcal{\ell }}\mathit{cos}\left(\varpi {\tau }_0\right){\varphi }_1-{\stackrel{\sim }{\alpha }}_2\frac{\partial {\varphi }_1}{\partial {\tau }_0}+\left(-z_1{\varphi }_2-z_2{\varphi }_1+\frac{{{\varphi }_1}^3}{6}\right){\omega }^2-2\left(\frac{\partial z_1}{\partial {\tau }_1}+\frac{\partial z_2}{\partial {\tau }_0}+z_1\frac{\partial z_1}{\partial {\tau }_0}\right)\frac{\partial {\varphi }_1}{\partial {\tau }_0}-\frac{{\partial }^2{\varphi }_1}{\partial {{\tau }_1}^2}-\frac{{\partial }^2{\varphi }_1}{\partial {{\tau }_0}^2}\left({z_1}^2+2z_2\right)-2\frac{\partial z_1}{\partial {\tau }_0}\left(\frac{\partial {\varphi }_2}{\partial {\tau }_0}+\frac{\partial {\varphi }_1}{\partial {\tau }_1}\right)-2z_1\left(\frac{{\partial }^2{\varphi }_2}{\partial {{\tau }_0}^2}+2\frac{{\partial }^2{\varphi }_1}{\partial {\tau }_0\partial {\tau }_1}\right)-2\left(\frac{{\partial }^2{\varphi }_2}{\partial {\tau }_0\partial {\tau }_1}+\frac{{\partial }^2{\varphi }_1}{\partial {\tau }_0\partial {\tau }_2}\right).$

(13)

The previous equations consist of partial differential equations (PDEs) in three groups. These equations are from second order that are solvable sequentially which promote the significance of the solutions of the first group (11). Consequently, they own the below forms in their general solutions

$z_1=A_1{e}^{i{\tau }_0}+{\overline{A}}_1{e}^{-i{\tau }_0},{\varphi }_1=A_2{e}^{i\omega {\tau }_0}+{\overline{A}}_2{e}^{-i\omega {\tau }_0};i=\sqrt{-1}.$

(14)

where $A_1$ and $A_2$ are complex functions that can later be determined, and they are given as functions of slow time scales ${\tau }_1$ and ${\tau }_2$ .

Substituting (14) into (12), and then omitting the terms that lead to secular terms to avoid appearance. Therefore, we can get

$\frac{\partial A_1}{\partial {\tau }_1}=0,\frac{\partial A_2}{\partial {\tau }_1}=0.$

(15)

Using both of the solutions (14) and the last conditions (15), the second order’s solutions can be obtained as follows

$z_2={\omega }^2{A}_2{\overline{A}}_2+\frac{{\varpi }^2\stackrel{\sim }{a}}{2\mathcal{\ell }\left(1-{\varpi }^2\right)}{e}^{i\varpi {\tau }_0}-\frac{3{\omega }^2{A_2}^2}{2\left(1-4{\omega }^2\right)}{e}^{2i\omega {\tau }_0}+CC,{\varphi }_2=\frac{-i{\varpi }^2\stackrel{\sim }{b}}{2\mathcal{\ell }\left({\omega }^2-{\varpi }^2\right)}{e}^{i\varpi {\tau }_0}-\frac{\omega \left(\omega +2\right)A_1{A}_2}{2\omega +1}{e}^{i\left(\omega +1\right){\tau }_0}+\frac{\omega \left(\omega -2\right){\overline{A}}_1{A}_2}{2\omega -1}{e}^{i\left(\omega -1\right){\tau }_0}+CC.$

(16)

here $CC$ indicates for the preceding terms’ complex conjugate.

Inserting (14)–(16) in equation (13) of the third order and eliminating the secular terms to acquire the following conditions

$i{\stackrel{\sim }{\alpha }}_1{A}_1+3\stackrel{\sim }{\beta }{y_r}^2{A}_1+2i\frac{\partial A_1}{\partial {\tau }_2}+\frac{6{\omega }^2{A}_1{A}_2{\overline{A}}_2}{\left(4{\omega }^2-1\right)}\left({\omega }^2-1\right)=0,i\omega {\stackrel{\sim }{\alpha }}_2{A}_2+2i\omega \frac{\partial A_2}{\partial {\tau }_2}+\frac{6{\omega }^2\left({\omega }^2-1\right)}{\left(4{\omega }^2-1\right)}{A}_1{\overline{A}}_1{A}_2+\frac{{\omega }^2\left(1+7{\omega }^2-8{\omega }^4\right)}{2\left(4{\omega }^2-1\right)}{A_2}^2{\overline{A}}_2=0.$

(17)

thus, the third-order solutions possess the following forms

$z_3=\frac{i\stackrel{\sim }{b}{\varpi }^3\left(\varpi +2\omega \right)A_2}{2\mathcal{\ell }\left({\omega }^2-{\varpi }^2\right)\left(1-{\left(\omega +\varpi \right)}^2\right)}{e}^{i\left(\omega +\varpi \right){\tau }_0}-\frac{3\omega \left(\omega +1\right)A_1{A_2}^2}{4\left(2\omega +1\right)}{e}^{i\left(2\omega +1\right){\tau }_0}+\frac{i\stackrel{\sim }{b}{\varpi }^3\left(\varpi -2\omega \right){\overline{A}}_2}{2\mathcal{\ell }\left({\omega }^2-{\varpi }^2\right)\left(1-{\left(\varpi -\omega \right)}^2\right)}{e}^{i\left(\varpi -\omega \right){\tau }_0}+\frac{3\omega \left(\omega -1\right){\overline{A}}_1{A_2}^2}{4\left(2\omega -1\right)}{e}^{i\left(2\omega -1\right){\tau }_0}+CC,{\varphi }_3=\frac{{\stackrel{\sim }{f}}_1}{2\left({\omega }^2-{p_1}^2\right)}{e}^{i{p}_1{\tau }_0}-\frac{i\stackrel{\sim }{b}{\varpi }^3\left(\varpi +2\right)A_1}{2\mathcal{\ell }\left({\omega }^2-{\varpi }^2\right)\left({\omega }^2-{\left(\varpi +1\right)}^2\right)}{e}^{i\left(\varpi +1\right){\tau }_0}-\frac{i\stackrel{\sim }{b}{\varpi }^3\left(\varpi -2\right){\overline{A}}_1}{2\mathcal{\ell }\left({\omega }^2-{\varpi }^2\right)\left({\omega }^2-{\left(\varpi -1\right)}^2\right)}{e}^{i\left(\varpi -1\right){\tau }_0}-\frac{\stackrel{\sim }{a}\varpi \left[{\left(\varpi +\omega \right)}^2-1\right]A_2}{2\mathcal{\ell }\left(1-{\varpi }^2\right)\left(2\omega +\varpi \right)}{e}^{i\left(\varpi +\omega \right){\tau }_0}+\frac{\stackrel{\sim }{a}\varpi \left[{\left(\varpi -\omega \right)}^2-1\right]{\overline{A}}_2}{2\mathcal{\ell }\left(1-{\varpi }^2\right)\left(2\omega -\varpi \right)}{e}^{i\left(\varpi -\omega \right){\tau }_0}+\frac{\omega \left(\omega +1\right)\left(\omega +2\right)\left(\omega +3\right){A_1}^2{A}_2}{4\left(2\omega +1\right)\left(\omega +1\right)}{e}^{i\left(\omega +2\right){\tau }_0}+\frac{\omega \left(\omega -1\right)\left(\omega -2\right)\left(\omega -3\right){{\overline{A}}_1}^2{A}_2}{4\left(2\omega -1\right)\left(\omega -1\right)}{e}^{i\left(\omega -2\right){\tau }_0}+\frac{{\omega }^2\left(49{\omega }^2-1\right){A_2}^3}{6\left(4{\omega }^2-1\right)\left({\omega }^2-9\right)}{e}^{3i{\tau }_0}+CC.$

(18)

the unknown functions $A_1$ and $A_2$ can be specified according to the conditions of deleting secular terms (15) and (17) besides the beginning prerequisites listed below

$\phi \left(0\right)={\eta }_{01},y\left(0\right)={\eta }_{03},\dot{\phi }\left(0\right)={\eta }_{02},\dot{y}\left(0\right)={\eta }_{04},$

(19)

where ${\eta }_{01},{\eta }_{02},{\eta }_{03},$ and ${\eta }_{04}$ are known quantities.

It must be noted that the obtained analytic solutions are fundamental to enhancing predictability, control, and performance while ensuring safety and efficiency. They provide a robust framework for understanding and managing the behavior of complex systems, leading to advancements across multiple fields and technologies. Here are some key benefits:

Knowing the solutions allows engineers to design systems that operate within desired parameters, avoiding resonant frequencies that could lead to excessive vibrations and potential failure. The obtained outcomes enable the development of control strategies to maintain desired motion patterns, enhancing system performance and safety. Stability analysis of these solutions helps in identifying safe operating conditions, ensuring that the system remains within stable regions under varying loads and disturbances. In mechanical systems, operating at such solutions can minimize energy consumption by reducing unnecessary or chaotic vibrations. In applications like automotive and aerospace engineering, understanding and controlling solutions help reduce noise generated by mechanical vibrations, leading to quieter operation.

Conditions for solvability and modulation equations

The obtained analytic solutions are acceptable when their denominators have non-zero value, that is, when the vibrations escape from resonance cases. ²⁴ These situations arise whenever one or more of the denominators in the second and third approximations equals zero. So, the detected resonance can be selected as follows:

• When $p_1=\omega$ is met, primary (main) external resonance takes place.

• When $\omega =\varpi ,\varpi =1,\omega =1,2\omega =1,2\omega =\varpi ,$ and $\omega =3$ are realized, the case of internal resonance is presented.

• At $\omega +\varpi =1,\omega =\varpi \pm 1$ , combined resonances are obtained.

If one of the aforementioned resonance situations is satisfied, the behavior of the system under investigation will be highly challenging. Then, we must modify the technique being used.

We’ll look into primary external and internal resonances that happen simultaneously in order to address this state. As a result, we consider the occurrence of the following two situations at once:

$\omega \approx 1,p_1\approx \omega .$

(20)

These relations show how near $\omega$ and $p_1$ are to $1$ and $\omega$ . To fulfill this purpose, introducing dimensionless amounts is required. They are known as detuning parameters ${\sigma }_1$ and ${\sigma }_2$ (which specify how far vibrations are from the strict resonance) in the way

$\omega =1+{\sigma }_1,p_1=\omega +{\sigma }_2.$

(21)

Then, we may assume they are of $O\left(\epsilon \right)$ as following

${\sigma }_j=\epsilon {\stackrel{\sim }{\sigma }}_j,j=1,2.$

(22)

Let us substitute (21) and (22) into (12) and (13), then eliminate terms that give secular ones to have the conditions of solvability for the second and third orders as follows

$\frac{\partial A_1}{\partial {\tau }_1}=0,\frac{\partial A_2}{\partial {\tau }_1}=0,i{\stackrel{\sim }{\alpha }}_1{A}_1+3\stackrel{\sim }{\beta }{y_r}^2{A}_1+2i\frac{\partial A_1}{\partial {\tau }_2}+\frac{6{\omega }^2\left({\omega }^2-1\right)A_1{A}_2{\overline{A}}_2}{4{\omega }^2-1}+\frac{3{\omega }^2{\left(\omega -1\right)}^2{\overline{A}}_1{A_2}^2}{2\omega -1}{e}^{2i{\stackrel{\sim }{\sigma }}_1{\tau }_1}=0,\frac{{\stackrel{\sim }{f}}_1}{2}{e}^{i{\stackrel{\sim }{\sigma }}_2{\tau }_1}-i\omega {\stackrel{\sim }{\alpha }}_2{A}_2-2i\omega \frac{\partial A_2}{\partial {\tau }_2}-\frac{6{\omega }^2\left({\omega }^2-1\right)}{\left(4{\omega }^2-1\right)}{A}_1{\overline{A}}_1{A}_2-\frac{{\omega }^2\left(1+7{\omega }^2-8{\omega }^4\right)}{2\left(4{\omega }^2-1\right)}{A_2}^2{\overline{A}}_2=0.$

(23)

The preceding solvability conditions distinguish that they are generating a system of four nonlinear PDE regarding with the undetermined functions $A_1$ and $A_2$ . Furthermore, the functions are depending upon the scale ${\tau }_2$ only, as obvious from condition (23). For this reason, these functions are presented in a polar form as follows

$A_j=\frac{1}{2}{\stackrel{\sim }{h}}_j\left({\tau }_2\right)e^{i{\psi }_j\left({\tau }_2\right)},h_j=\epsilon {\stackrel{\sim }{h}}_j\left(j=1,2\right),$

(24)

where ${\psi }_j$ and ${\stackrel{\sim }{h}}_j$ express the phases and amplitudes of the solutions $z$ and $\varphi$ , respectively, in which they are real functions.

The simulations above allow us to rewrite derivative operators of the first order for $A_1$ and $A_2$ in the following ways

$\frac{d{A}_j}{d\tau }={\epsilon }^2\frac{\partial A_j}{\partial {\tau }_2},$

(25)

The PDE (23) can be transformed into ordinary differential equations (ODEs) using the formula (25). Inserting (24), (25), (8), and (10) in addition to the next adjusted phases

${\theta }_1\left({\tau }_2\right)=2\left({\stackrel{\sim }{\sigma }}_1{\tau }_1+{\psi }_2\left({\tau }_2\right)-{\psi }_1\left({\tau }_2\right)\right),{\theta }_2\left({\tau }_2\right)={\stackrel{\sim }{\sigma }}_2{\tau }_1-{\psi }_2\left({\tau }_2\right),$

(26)

into (23), by detaching the real and imaginary parts. It is easy to get the below autonomous system of ODE regarding to $h_1,h_2,{\theta }_1,$ and ${\theta }_2$

$\frac{d{h}_1}{d\tau }=-\frac{1}{2}{\alpha }_1{h}_1-\frac{3{\omega }^2{\left(\omega -1\right)}^2{h}_1{h_2}^2}{8\left(2\omega -1\right)}\mathit{sin}\hspace{0.17em}{\theta }_1,\frac{d{\theta }_1}{d\tau }=2\left({\sigma }_1+{\sigma }_2-\frac{d{\theta }_2}{d\tau }\right)-3\beta y_r^2-\frac{3{\omega }^2\left({\omega }^2-1\right)}{2\left(4{\omega }^2-1\right)}{h}_2^2-\frac{3{\omega }^2{\left(\omega -1\right)}^2{h_2}^2}{4\left(2\omega -1\right)}\mathit{cos}\hspace{0.17em}{\theta }_1,\frac{d{h}_2}{d\tau }=\frac{f_1}{2\omega }\mathit{sin}\hspace{0.17em}{\theta }_2-\frac{{\alpha }_2}{2}{h}_2,\frac{d{\theta }_2}{d\tau }=\frac{f_1}{2\omega h_2}\mathit{cos}\hspace{0.17em}{\theta }_2+{\sigma }_2-\frac{3\omega \left({\omega }^2-1\right)}{4\left(4{\omega }^2-1\right)}{h}_1^2-\frac{\omega \left(1+7{\omega }^2-8{\omega }^4\right)}{16\left(4{\omega }^2-1\right)}{h}_2^2.$

(27)

The above equations are known by the equations of modulation for the amplitudes $h_j$ and adjusted phases ${\theta }_j\left(j=1,2\right)$ , where it is made up of four first-order ODEs that can be gained when two resonance cases happen at once. Based on the below data, these equations are solved numerically, to show their behaviors at certain values of the parameters of the system, as seen in Figures 2–7

$F_0=0.1N,c_1=0.30N.m.s/rad,c_2=0.20N.s/m,k_1=20kg.s^{-2},k_2=70kg.s^{-2},a=0.3m,b=0.5m,\Omega =2.2rad.s^{-1},{\Omega }_1=2.5rad.s^{-1},{\sigma }_1=\omega -1,{\sigma }_2=p_1-\omega ,\mathcal{\ell }=1m,m_1=2kg,\epsilon =0.005,{\omega }_1\left(=3.65,3.16,2.82\right)s^{-1},{\omega }_2\left(=4.42,3.74,3.13\right)s^{-1}.$

to examine the importance of the used parameters on the behavior of dynamical systems, Figures 2–7 are sketched to present the transient responses of the adjusted amplitudes $h_j$ and phases ${\theta }_j\left(j=1,2\right)$ at specified values of the frequencies ${\omega }_j\left(j=1,2\right)$ . These values line up with the pendulum arm’s typical length $\mathcal{\ell }$ , the damping coefficients $c_j\left(j=1,2\right)$ , and the mass $m_1$ , as graphed in portions of these figures. Curves of them are decaying with $\tau$ , which express their stability behavior. Figure 2 plots the variation of $h_1$ with time when ${\omega }_j$ and $c_j$ have the values $\left(3.65,3.16,2.82\right),\left(4.42,3.74,3.13\right),$ and $\left(0.03,0.15,0.30\right),$ $\left(0.02,0.10,0.20\right)$ , respectively. According to the mentioned different values of parameters, the variation is presented in Figure 2(a)–(d). It is obvious that the drawn curves in these figures decrease over time which is in an agreement with the first equation (27).

Figure 2.

Variation of $h_1$ versus $\tau$ at (a) ${\mathit{\omega }}_{\mathbf{1}}\left(=\mathbf{3.65},\mathbf{3.16},\mathbf{2.82}\right),$ , (b) ${\mathit{\omega }}_{\mathbf{2}}\left(=\mathbf{4.42},\mathbf{3.74},\mathbf{3.13}\right),\left(\mathit{c}\right)\hspace{0.17em}{\mathit{c}}_{\mathbf{1}}\left(=\mathbf{0.03},\mathbf{0.15},\mathbf{0.30}\right),$ and (d) ${\mathit{c}}_{\mathbf{2}}\left(=\mathbf{0.02},\mathbf{0.10},\mathbf{0.20}\right)$ .

Figure 3.

The time histories of $h_2$ at distinct values of ${\omega }_j\left(j=1,2\right)$ and $c_j$ .

Figure 4.

Phase ${\mathit{\theta }}_{\mathbf{1}}\left(\mathit{\tau }\right)$ (a) when ${\mathit{\omega }}_{\mathbf{1}}\left(=\mathbf{3.65},\mathbf{3.16},\mathbf{2.82}\right),$ (b) when ${\mathit{\omega }}_{\mathbf{2}}\left(=\mathbf{4.42},\mathbf{3.74},\mathbf{3.13}\right),$ (c) when ${\mathit{c}}_{\mathbf{1}}\left(=\mathbf{0.03},\mathbf{0.15},\mathbf{0.30}\right),$ and (d) when ${\mathit{c}}_{\mathbf{2}}\left(=\mathbf{0.02},\mathbf{0.10},\mathbf{0.20}\right)$ .

Figure 5.

Displays ${\mathit{\theta }}_{\mathbf{2}}\left(\mathit{\tau }\right)$ for various values of ${\mathit{\omega }}_{\mathit{j}}$ and ${\mathit{c}}_{\mathit{j}}$ .

Figure 6.

Curves in the plane ${\theta }_1{h}_1$ at selected values of ${\omega }_j$ and $c_j$ .

Figure 7.

Phase plane curves ${\theta }_2{h}_2$ for distinct values of ${\omega }_j$ and $c_j$ at (a) ${\mathit{\omega }}_{\mathbf{1}}\left(=\mathbf{3.65},\mathbf{3.16},\mathbf{2.82}\right),$ (b) ${\mathit{\omega }}_{\mathbf{2}}\left(=\mathbf{4.42},\mathbf{3.74},\mathbf{3.13}\right),$ (c) ${\mathit{c}}_{\mathbf{1}}\left(=\mathbf{0.03},\mathbf{0.15},\mathbf{0.30}\right),$ and (d) ${\mathit{c}}_{\mathbf{2}}\left(=\mathbf{0.02},\mathbf{0.10},\mathbf{0.20}\right)$ .

Examining the indicated curves in Figure (3) reveals that when ${\omega }_j$ and $c_j$ have the various selected values mentioned above, the amplitude $h_2$ changes with $\tau$ . Generally speaking, the values of these parameters affect the way that all drawn waves decay, and we can note that if ${\omega }_1$ increases, the amplitude of wave increases and the number of oscillations decreases, on the contrary, the performance of ${\omega }_2$ . The plotted curve in portion (c) shows that the change of the damping coefficient $c_1$ doesn’t affect the behavior of the amplitude $h_2$ . We may see in portion (d) of the same figure that when the coefficient of damping $c_2$ increases, the amplitude of $h_2$ decreases while the oscillations’ number remains steady.

On the reverse side, Figures 4 and 5 show the time change of the adjusted phases ${\theta }_j\left(\tau \right)$ when ${\omega }_j$ and $c_j$ have distinct values. It is observed that ${\theta }_1$ value changes have a negligible effect, as graphed in Figure 4(a), (c), and (d), while this alteration is progressively obvious with various values of ${\omega }_2$ , as indicated in Figure 4(b). The behavior of ${\theta }_2$ waves has a positive impact with change of the values of ${\omega }_1,{\omega }_2,$ and $c_2$ as explored, respectively, in the parts (a), (b), and (d) of Figure 5. However, ${\theta }_2$ has no variation, to some extent, with the values of $c_1$ , as explored in Figure 5(c). The reason is owing to the first and fourth equations of system (27), in which the last equation doesn’t depend explicitly on $c_1$ , while it depends upon $c_2$ directly, as seen in the dimensionless parameters (5).

Figures 6 and 7 display the phase plane charts of the amplitudes ${\mathit{h}}_{\mathit{j}}\left(\mathit{j}=\mathbf{1},\mathbf{2}\right)$ and the adjusted phases ${\mathit{\theta }}_{\mathit{j}}$ for different values of ${\omega }_j$ and $c_j$ . It is significant to note that the drawn curves in the plane ${\mathit{\theta }}_{\mathbf{1}}{\mathit{h}}_{\mathbf{1}}$ have a decaying behavior with time, in which ${\mathit{\omega }}_{\mathbf{2}}$ and ${\mathit{c}}_{\mathbf{1}}$ have a positive impact on these curves, as graphed in Figure 6(b) and (c). Curves in Figure 6(a) and (d) don’t vary with the changes of ${\mathit{\omega }}_{\mathbf{1}}$ and ${\mathit{c}}_{\mathbf{2}}$ . The inspection of the curves of Figure 7 reveals spiral arcs pointing in one direction, which demonstrates the stability of these amplitudes and phases. The sensitivity of the dynamical system (17) to the chosen values ${\mathit{\omega }}_{\mathbf{1}},{\mathit{\omega }}_{\mathbf{2}},$ and ${\mathit{c}}_{\mathbf{2}}$ is much higher than the value of ${\mathit{c}}_{\mathbf{1}}$ , as seen from the plotted curves in Figure 7.

In Figures 8 and 9, the approximate solutions $\phi$ and $y$ through time $\tau$ have been drawn to represent the time-historical behavior pertaining to these solutions. These figures were computed using earlier data.

Figure 8.

Temporal time behavior of $\phi$ when ${\omega }_j$ and $c_j$ vary.

Figure 9.

Impact of changing ${\mathit{\omega }}_{\mathit{j}}$ and ${\mathit{c}}_{\mathit{j}}$ on the solution $\mathit{y}$ .

In Figure 8, the variation of the rotation angle $\mathit{\phi }$ is described as time $\mathit{\tau }$ passes, while Figure 9 displays the variation of spring’s extension $\mathit{y}$ with time $\mathit{\tau }$ when the parameters being utilized are the same values. Look at the graphs in portions of Figure 8, where it is evident that as ${\mathit{\omega }}_{\mathbf{1}}$ rises, the plotted waves having periodic forms and their amplitudes decrease, as seen in portion (a), while portion (b) reveals decrease of the amplitude, to some extent, when ${\mathit{\omega }}_{\mathbf{2}}$ increases. On the other hand, Figure 9(a) demonstrates a reduction in wave amplitudes when ${\mathit{\omega }}_{\mathbf{1}}$ decreases, whereas an increment of the wave amplitudes is found in Figure 9(b) with the raise of ${\mathit{\omega }}_{\mathbf{2}}$ values. The change of the approximate solutions with the values of ${\mathit{c}}_{\mathbf{1}}$ and ${\mathit{c}}_{\mathbf{2}}$ seems slightly, as represented in portions (c) and (d) of Figures 8 and 9. According to these curves, the solutions are considered to be stable

Finally, Figures 10 and 11 display the trajectory projection of the modulation system on the phase plane $\phi {\phi }^{\prime }$ and $y{y}^{\prime }$ . These waves are obtained at different parameters ${\mathit{\omega }}_{\mathit{j}}$ and $c_j$ , respectively, in which their curves have closed forms, that are consistent with the aforementioned discussion and assert the stability of the obtained solutions.

Figure 10.

Curves’ phase plane of $\phi$ at the considered values of ${\omega }_j$ and $c_j$ .

Figure 11.

Curves’ phase plane of $y$ at the considered values of ${\omega }_j$ and $c_j$ .

The algorithm of fourth-order Runge–Kutta is applied to gain the numerical solutions (NSs) of the basic EOM. They are graphed at $c_1=0.30,c_2=0.20,{\omega }_1=3.16,{\omega }_2=3.74$ in Figure 12(a) and (b) for the solutions $\phi$ and $y$ , respectively. On the other hand, the corresponding analytical solutions (ASs), are drawn for the aforementioned parameters in Figure 12(c) and (d). Figure 12(e) and (f) show a comparison between the NS and AS that demonstrates their great compatibility which demonstrates that the utilized perturbation approach is accurate.

Figure 12.

Comparison between the AS and NS at ${\mathit{\omega }}_{\mathbf{1}}=\mathbf{3.16}$ and ${\mathit{\omega }}_{\mathbf{2}}=\mathbf{3.74}$ .

Steady-state solutions

The purpose of this section is to look the oscillations of the examined dynamical system at the scenario of steady-state. In general, steady-state oscillations will be produced when transient processes terminate under the impact of dampening. The amplitudes and adjusted phases of steady-state may be attained in light of zero values of the system (27), that is, at $\frac{d{h}_j}{d\tau }=\frac{d{\theta }_j}{d\tau }=0\left(j=1,2\right).$ Thus, we are able to write

${\alpha }_1{h}_1+\frac{3{\omega }^2{\left(\omega -1\right)}^2{h}_1{h_2}^2}{4\left(2\omega -1\right)}\mathit{sin}\hspace{0.17em}{\theta }_1=0,2\left({\sigma }_1+{\sigma }_2\right)-3\beta y_r^2-\frac{3{\omega }^2\left({\omega }^2-1\right)}{2\left(4{\omega }^2-1\right)}{h}_2^2-\frac{3{\omega }^2{\left(\omega -1\right)}^2{h_2}^2}{4\left(2\omega -1\right)}\mathit{cos}\hspace{0.17em}{\theta }_1=0,f_1\hspace{0.17em}\mathit{sin}\hspace{0.17em}{\theta }_2-\omega {\alpha }_2{h}_2=0,f_1\hspace{0.17em}\mathit{cos}\hspace{0.17em}{\theta }_2+2\omega h_2{\sigma }_2-\frac{3{\omega }^2\left({\omega }^2-1\right)}{2\left(4{\omega }^2-1\right)}{h}_1^2{h}_2-\frac{{\omega }^2\left(1+7{\omega }^2-8{\omega }^4\right)}{8\left(4{\omega }^2-1\right)}{h}_2^3=0.$

(28)

This system is composed of four algebraic equations that are expressed in terms of the amplitudes $h_j$ and adjusted phases ${\theta }_j$ . When the modified phases ${\theta }_j$ are considered to be absent, system (28) generates two nonlinear algebraic equations between the amplitudes $h_j$ and the parameters of detuning ${\sigma }_j$ , as follows

${\left[\frac{3{\omega }^2{\left(\omega -1\right)}^2{h}_1{h_2}^2}{4\left(2\omega -1\right)}\right]}^2={{\alpha }_1}^2{h_1}^2+{\left[2h_1\left({\sigma }_1+{\sigma }_2\right)-3h_1\beta y_r^2-\frac{3{\omega }^2\left({\omega }^2-1\right)}{2\left(4{\omega }^2-1\right)}{h}_1{h}_2^2\right]}^2,{f_1}^2={\omega }^2{{\alpha }_2}^2{h_2}^2+4{\omega }^2{\left[\frac{3\omega \left({\omega }^2-1\right)}{4\left(4{\omega }^2-1\right)}{h}_1^2{h}_2+\frac{\omega \left(1+7{\omega }^2-8{\omega }^4\right)}{16\left(4{\omega }^2-1\right)}{h}_2^3-{\sigma }_2{h}_2\right]}^2.$

(29)

It is important to note that $h_1$ and $h_2$ represent, respectively, the amplitudes of the longitudinal and swing oscillations.

In light of the equations of system (29), the stability of the FP can be assessed. Therefore, curves of Figure 13 are plotted in the plane $h_2{h}_1$ at $F_0=2,k_1=20,$ and ${\omega }_2=2.47$ . The intersection of these curves generates FP, in which they are zoomed and marked by small gray circles, as in Figure 14. Figure 14(a) and (c) contain one unstable fixed point, which doesn’t satisfy all the conditions of (34), while three unstable FPs are indicated in Figure 13(b).

Figure 13.

Represented curves of the system (29).

Figure 14.

Unstable FP in the plane $h_2{h}_1$ .

To accomplish the system’s performance in a region that is quite near the FP, consequently, the following expressions of $h_j$ and ${\theta }_j$ are taken into consideration in the system (27)

$h_1=h_{10}+h_{11},{\theta }_1={\theta }_{10}+{\theta }_{11},h_2=h_{20}+h_{21},{\theta }_2={\theta }_{20}+{\theta }_{21},$

(30)

here $h_{j0}$ and ${\theta }_{j0}$ represent the solutions at the scenario of steady-state, while $h_{j1}$ and ${\theta }_{j1}$ denote the very small perturbations regarding $h_{j0}$ and ${\theta }_{j0}$ , respectively. The substitution of (30) in (27) yields the following form

$\frac{d{h}_{11}}{d\tau }=-\frac{{\alpha }_1}{2}{h}_{11}-\frac{3{\omega }^2{\left(\omega -1\right)}^2}{8\left(2\omega -1\right)}[h_{10}{h_{20}}^2\hspace{0.17em}\mathit{cos}\left({\theta }_{10}\right){\theta }_{11}+2h_{10}{h}_{20}{h}_{21}\hspace{0.17em}\mathit{sin}\left({\theta }_{10}\right)+h_{11}{h_{20}}^2\hspace{0.17em}\mathit{sin}\left({\theta }_{10}\right)],h_{10}\frac{d{\theta }_{11}}{d\tau }+2h_{10}\frac{d{\theta }_{21}}{d\tau }=-3\beta {y_r}^2{h}_{11}+2\left({\sigma }_1+{\sigma }_2\right)h_{11}-\frac{3{\omega }^2\left({\omega }^2-1\right)}{2\left(4{\omega }^2-1\right)}\left[2h_{10}{h}_{20}{h}_{21}+h_{11}{h_{20}}^2\right]-\frac{3{\omega }^2{\left(\omega -1\right)}^2}{4\left(2\omega -1\right)}\left[-h_{10}{h_{20}}^2\hspace{0.17em}\mathit{sin}\left({\theta }_{10}\right){\theta }_{11}+2h_{10}{h}_{20}{h}_{21}\hspace{0.17em}\mathit{cos}\left({\theta }_{10}\right)+h_{11}{h_{20}}^2\hspace{0.17em}\mathit{cos}\left({\theta }_{10}\right)\right],\frac{d{h}_{21}}{d\tau }=-\frac{{\alpha }_2}{2}{h}_{21}+\frac{f_1}{2\omega }\mathit{cos}\left({\theta }_{20}\right){\theta }_{21},h_{20}\frac{d{\theta }_{21}}{d\tau }={\sigma }_2{h}_{21}-\frac{3\omega \left({\omega }^2-1\right)}{4\left(4{\omega }^2-1\right)}\left[{h_{10}}^2{h}_{21}+2h_{10}{h}_{20}{h}_{11}\right]-\frac{3\omega \left(1+7{\omega }^2-8{\omega }^4\right)}{16\left(4{\omega }^2-1\right)}{h_{20}}^2{h}_{21}-\frac{f_1}{2\omega }\mathit{sin}\left({\theta }_{20}\right){\theta }_{21}.$

(31)

Based on the unknown perturbed functions ${\theta }_{j1}$ and $h_{j1}$ of the adjusted phases and amplitudes, their solutions can be represented as linear superposition of $k_s{e}^{\lambda \tau }\left(s=1,2,3,4\right)$ where $k_s$ and $\lambda$ are eigenvalues of the constants and unknown perturbation, respectively. If the solutions ${\theta }_{j0}$ and $h_{j0}$ are stable asymptotically, then the real parts of the roots of the below characteristic equation of (31) should be negative

${\lambda }^4+{\Gamma }_1{\lambda }^3+{\Gamma }_2{\lambda }^2+{\Gamma }_3\lambda +{\Gamma }_4=0,$

(32)

where

${\Gamma }_1=\frac{1}{2}\left[\frac{f_1\hspace{0.17em}\mathit{sin}\hspace{0.17em}{\theta }_{20}}{\omega h_{20}}+d_0{h}_{20}^2\hspace{0.17em}\mathit{sin}\hspace{0.17em}{\theta }_{10}+{\alpha }_1+{\alpha }_2\right],{\Gamma }_2=\frac{1}{4}\{-2d_0^2{h}_{20}^4+{\alpha }_1{\alpha }_2+d_0\hspace{0.17em}\mathit{sin}\hspace{0.17em}{\theta }_{10}{h}_{20}^2\left({\alpha }_2-{\alpha }_1\right)+\frac{f_1}{\omega h_{20}}[2d_2\hspace{0.17em}\mathit{cos}\hspace{0.17em}{\theta }_{20}{h}_{10}^2+\left(6d_3\hspace{0.17em}\mathit{cos}\hspace{0.17em}{\theta }_{20}+d_0\hspace{0.17em}\mathit{sin}\hspace{0.17em}{\theta }_{10}\hspace{0.17em}\mathit{sin}\hspace{0.17em}{\theta }_{20}\right)h_{20}^2+\mathit{sin}\hspace{0.17em}{\theta }_{20}\left({\alpha }_1+{\alpha }_2\right)-2\hspace{0.17em}\mathit{cos}\hspace{0.17em}{\theta }_{20}{\sigma }_2]+4d_0\hspace{0.17em}\mathit{cos}\hspace{0.17em}{\theta }_{10}{h}_{20}^2\left[4d_2{h}_{10}^2-d_1{h}_{20}^2-3\beta y_r^2+2\left({\sigma }_1+{\sigma }_2\right)\right]\},$

${\Gamma }_3=\frac{1}{8h_{20}\omega }\{f_1[-2d_0^2{\mathit{sin}\hspace{0.17em}\theta }_{20}{h}_{20}^4+2d_2{\mathit{cos}\hspace{0.17em}\theta }_{20}{h}_{10}^2{\alpha }_1+6d_3{\mathit{cos}\hspace{0.17em}\theta }_{20}{h}_{20}^2{\alpha }_1+{\mathit{sin}\hspace{0.17em}\theta }_{20}{\alpha }_1{\alpha }_2-2{\mathit{cos}\hspace{0.17em}\theta }_{20}{\alpha }_1{\alpha }_2+d_0{\mathit{sin}\hspace{0.17em}\theta }_{10}{h}_{20}^2\left[{\mathit{sin}\hspace{0.17em}\theta }_{20}\left({\alpha }_2-{\alpha }_1\right)-2{\mathit{cos}\hspace{0.17em}\theta }_{20}\left(7d_2{h}_{10}^2-3d_3{h}_{20}^2+{\sigma }_2\right)\right]-4d_0{\mathit{cos}\hspace{0.17em}\theta }_{10}{\mathit{sin}\hspace{0.17em}\theta }_{20}{h}_{20}^2\left(d_1{h}_{20}^2+3\beta y_r^2-2\left({\sigma }_1+{\sigma }_2\right)\right)]+d_0\omega h_{20}^3{\alpha }_2(-2d_0{h}_{20}^2-{\mathit{sin}\hspace{0.17em}\theta }_{10}{\alpha }_1+4{\mathit{cos}\hspace{0.17em}\theta }_{10}\left(4d_2{h}_{10}^2-d_1{h}_{20}^2-3\beta y_r^2+2\left({\sigma }_1+{\sigma }_2\right)\right))\},$

${\Gamma }_4=\frac{d_0{h}_{20}{f}_1}{16\omega }\{-\left[\mathit{sin}\hspace{0.17em}{\theta }_{20}{\alpha }_2+2\hspace{0.17em}\mathit{cos}\hspace{0.17em}{\theta }_{20}\left(3d_3{h}_{20}^2-{\sigma }_2\right)\right](2d_0{h}_{20}^2+\mathit{sin}\hspace{0.17em}{\theta }_{10}{\alpha }_1++4\hspace{0.17em}\mathit{cos}\hspace{0.17em}{\theta }_{10}\times \left[d_1{h}_{20}^2+3\beta y_r^2-2\left({\sigma }_1+{\sigma }_2\right)\right])+2d_2{\mathit{cos}\hspace{0.17em}\theta }_{20}{h}_{10}^2(2(d_0+\left(4d_0+6d_1\right)\mathit{cos}\hspace{0.17em}{\theta }_{10}-2d_0\hspace{0.17em}\mathit{cos}\hspace{0.17em}2{\theta }_{10}{h}_{20}^2-\hspace{0.17em}\mathit{sin}\hspace{0.17em}{\theta }_{10}{\alpha }_1+4\hspace{0.17em}\mathit{cos}\hspace{0.17em}{\theta }_{10}\left[2\left({\sigma }_1+{\sigma }_2\right)-3\beta y_r^2\right])\}.$

(33)

where

$d_0=\frac{3{\omega }^2{\left(\omega -1\right)}^2}{8\left(2\omega -1\right)},d_2=\frac{3\omega \left({\omega }^2-1\right)}{4\left(4{\omega }^2-1\right)},d_1=\frac{3{\omega }^2\left({\omega }^2-1\right)}{2\left(4{\omega }^2-1\right)},d_3=\frac{\omega \left(7{\omega }^2-8{\omega }^4+1\right)}{16\left(4{\omega }^2-1\right)}.$

The fundamental prerequisites for stability under the RHC ⁴⁸ for the specific steady state solutions are

${\Gamma }_1>0,{\Gamma }_1{\Gamma }_2-{\Gamma }_3>0,{\Gamma }_4\left({\Gamma }_1{\Gamma }_2{\Gamma }_3-{{\Gamma }_3}^2-{{\Gamma }_1}^2{\Gamma }_4\right)>0.{\Gamma }_1{\Gamma }_2{\Gamma }_3-{{\Gamma }_3}^2-{{\Gamma }_1}^2{\Gamma }_4>0,$

(34)

Assessment of the stability

The approach of non-linear stability is applied in the present section, to study the stability of the understudied system, and both the simulation and the application of the stability parameters are included. The frequency amplitudes ${\omega }_j\left(j=1,2\right)$ , detuning parameters ${\sigma }_j$ , and damping coefficient $c_j$ are among the most crucial factors in the stabilization process. Considering the previously selected values for these parameters, the stability diagrams of the system (31) are plotted in Figures 15–17. For various parametrical areas, the amplitudes $h_j$ are displayed, and the phase plane tracks in these areas display their properties. These figures show the frequency response curves with the varying values of ${\omega }_j$ and $c_j$ . The stability/instability areas for all possible FPs are examined through the range $-5\le {\sigma }_1\le 10$ of the curves of frequency responses. From Figure 15(a), it can be shown that, when ${\omega }_1=3.65$ , the stable fixed point is satisfied in the region $5.6<{\sigma }_1$ , while the domain of unstable FP is ${\sigma }_1\le 5.6$ . Similarly, the stable and unstable regions are obtained at ${\omega }_1\left(=3.16,2.82\right)$ are, respectively, $5.4<{\sigma }_1$ , $5.2<{\sigma }_1$ , ${\sigma }_1\le 5.4$ , and ${\sigma }_1\le 5.2$ . The dashed curves and the solid ones indicate the ranges of unstable and stable FPs, respectively. Each curve of this figure contains one critical point (the point between two distinct areas) and two peak points. When ${\omega }_2\left(=4.42,3.74,3.13\right)$ , the stability and instability areas have been graphed in the ranges $3.1<{\sigma }_2,$ $5.4<{\sigma }_2,$ $8.3<{\sigma }_2$ and ${\sigma }_2\le 3.1,$ ${\sigma }_2\le 5.4,$ ${\sigma }_2\le 8.3,$ respectively, as seen in Figure 15(b). One critical fixed point and two peak points may be found on each curve in this figure. It should be noticed from Figure 15(c) and (d) that the possible FPs lose their stabilities in the area ${\sigma }_1<5.2$ , while they are stable in the range $5.2\le {\sigma }_1$ .

Figure 15.

Stability/instability regions in ${\sigma }_1{h}_1$ at ${\sigma }_2=p_1-\omega$ when (a) ${\mathit{\omega }}_{\mathbf{1}}\left(=\mathbf{3.65},\mathbf{3.16},\mathbf{2.82}\right),\left(b\right)\hspace{0.17em}{\mathit{\omega }}_{\mathbf{2}}\left(=\mathbf{4.42},\mathbf{3.74},\mathbf{3.13}\right),\left(c\right)c_1\left(=0.03,0.15,0.30\right),{\mathit{c}}_{\mathbf{2}}\left(=\mathbf{0.02},\mathbf{0.10},\mathbf{0.20}\right)$ .

Figure 16.

Curves in the plane ${\sigma }_1{h}_2$ when ${\sigma }_2=p_1-\omega$ , at the same values of Figure 15.

Figure 17.

Undesired frequency response in ${\sigma }_1{h}_1$ and ${\sigma }_1{h}_2$ at ${\sigma }_2=0.0079$ when ${\mathit{\omega }}_{\mathbf{2}}=\mathbf{2.47}{\mathit{s}}^{-\mathbf{1}}$ and $F_0=2N$ .

Based on the previous analysis of Figure 15, we determine that the response curves in the plane ${\sigma }_1{h}_2$ have the same ranges of instability and stability FP, as drawn in portions of Figure 16. It is important to state that, when transcritical bifurcation occurs due to one fixed point, the connecting points between both curves characterize the wave’s critical values. From Figure 17(a) it is noticed that the area of unstable FP located at ${\omega }_1=2.47$ is ${\sigma }_1\le 4.9$ , while the stable one lies in the range $4.9<{\sigma }_1$ . One critical fixed point is observed, besides several peak ones on these curves, which it is called undesired frequency response.

Recognizing the stability and instability zones in vibrating dynamical systems yields a range of benefits across diverse applications. These include the following:

Recognizing stability zones enables engineers to design systems that function safely within these regions, reducing the likelihood of unforeseen failures. Understanding instability zones helps designers avoid parameters that cause excessive vibrations, resonance, and potential catastrophic failure. Maintaining operation within stability zones ensures smooth, efficient performance and minimizes mechanical wear. Awareness of these zones supports the tuning of system parameters for optimal performance and efficiency. Recognizing stability zones is crucial for designing effective vibration isolation systems that safeguard sensitive equipment from harmful vibrations. In structural engineering, stability analysis is applied to design buildings and infrastructure resilient to seismic activities by steering clear of instability zones that could amplify vibrations. Ensuring the stable operation of medical devices like heart pacemakers is vital for their dependable functioning.

Concluding remarks

The nonlinear planar motion of an influenced 2DOF damped elastic pendulum by harmonic force has been examined. The novelty of this work lies in the specific combination of damping, external forcing, and the elliptic path of the pendulum’s pivot point. The motion of this point is a unique aspect that distinguishes this study from other works on spring pendulums. The EOM have been obtained using Lagrangian equations. The novel analytic approximate solutions to these equations have been achieved utilizing the MSA to within three approximations. This approach helps capture the system’s behavior more accurately in the presence of nonlinearities. All resonance cases that have been detected are grouped, in which two of them are inspected simultaneously. This thorough examination of such cases adds depth to the understanding of the system’s dynamics. In line with the elimination of the secular terms, the solvability restrictions have been fulfilled, and the steady-state scenarios have been examined. These processes are considered crucial for understanding the long-term behavior of the system. The solutions' temporal trajectories, the resonance response curves, and the steady-state solutions have been displayed graphically. According to the RHC, all pertinent FPs at steady and unstable states have been located and shown. The stability zones have been looked and analyzed to determine the positive effects of various parameters, which provide practical insights for controlling and predicting the system’s responses.

Future Work

Future research could involve examining other dynamical models with 3DOF or 4DOF. More importantly, we could apply different perturbation methods to find analytical solutions. Additionally, we might model the rigid body problem as a pendulum, incorporating the movement of the support point along different paths.

Application

The examined model has gained prominence recently due to its industrial uses in seismic isolation systems for buildings and structures. In seismic engineering, a 2DOF vibrating pendulum system can be utilized as part of a seismic isolation system designed to protect buildings and infrastructure from earthquake-induced vibrations. The pendulum mechanism helps to absorb and dissipate seismic energy, reducing the amount of force transmitted to the structure. During an earthquake, the ground motion acts as an external harmonic force on the building. The distribution of mass and the structural layout can cause rotational moments that act on the building. The pendulum system can be tuned to counteract these moments, helping to stabilize the structure. The pendulum system allows for displacement in both horizontal and vertical directions, providing two degrees of freedom. This capability is essential for accommodating the complex and multi-directional nature of seismic waves.