The stability of 3-DOF triple-rigid-body pendulum system near resonances

In this article, the motion of three degree-of-freedom (DOF) dynamical system consisting of a triple rigid body pendulum (TRBP) in the presence of three harmonically external moments is studied. In view of the generalized coordinates of the system, Lagrange's equations are used to obtain the governing system of equations of motion (EOM). The analytic approximate solutions are gained up to the third approximation utilizing the approach of multiple scales (AMS) as novel solutions. The solvability conditions are determined in accordance with the elimination of secular terms. Therefore, the arising various resonances cases have been categorized and the equations of modulation have been achieved. The temporal histories of the obtained approximate solutions, as well as the resonance curves, are visually displayed to reveal the positive effects of the various parameters on the dynamical motion. The numerical results of the governing system are achieved using the fourth-order Runge–Kutta method. The visually depicted comparison of asymptotic and numerical solutions demonstrates high accuracy of the employed perturbation approach. The criteria of Routh–Hurwitz are used to investigate the stability and instability zones, which are then analyzed in terms of steady-state solutions. The strength of this work stems from its uses in engineering vibrational control applications which carry the investigated system a huge amount of importance.


Introduction
Due to their widespread applicability in everyday life, pendulum models have piqued the curiosity of many scholars in recent decades, particularly the first two decades of this century. However, the dynamics of many different types of pendulums were examined using various methods whether they are analytical or numerical.
In [1], the authors examined the instability positions of an inverted pendulum, in which if its pivot is harmonically pushed up and down with suitable frequency and amplitude, then its instability position can be abolished. The parametric pendulum is used in [2] to show how non-harmonic perturbations can be used to flip from one dynamic state to the other without changing the system's intrinsic characteristics. In [3], the authors investigated the motion of a simple pendulum as a dynamical model with single generalized coordinate to examine its periodic movement on an elliptic route. On the other hand, the chaotic motion of 2-DOF weakly nonlinear spring pendulum system for a fixed and moving suspension point is examined in [4] and [5], respectively, to reveal their chaotic behaviors near resonance with respect to the used parameters. The ASM [6] is used in [7] to explore the stability of a nonlinear spring pendulum, in which the improvement of this work is investigated in [8] when the pendulum is subjected to an external force.
Numerous publications, such as [9][10][11][12][13], have investigated various routes of the suspension points for the linear and nonlinear elastic pendulums with various DOFs. In [9], the suspension point of a damped elastic pendulum (EP) moves in a circular path. The obtained results are generalized in [10] for an elliptic route of this point, in which the approximate solutions (AS) of some special cases are obtained. In [11], the authors considered such a problem when the pivot point has a trajectory of Lissajous curve, while the general outcomes of this work are found in [12] for the planar motion of an EP connected with a rigid body to produce a dynamical model with 3 DOFs on the same route of its suspension point. A special case of these results is found in [13] for a fixed point of suspension. The cases of linear and nonlinear damped EP carrying a rigid body, when the pivot point has a path of ellipse, is examined in [14,15], respectively. The AMS is used in these works to gain the AS of the governing system of motion. Therefore, the requirements of solvability are obtained through the conditions of eliminating secular terms and the possible resonance cases are characterized. The vibrational motion of a rigid body regarding the equilibrium position is investigated numerically in [16] using the framework of ode45 solver of the Runge-Kutta method [17] from fourth order. Recently, a comparison between numerical solutions (NS) and the AS for the motions of rigid bodies pendulum is examined in [18,19] to highlight the good reliability between them and to explore the high accuracy of the adopted perturbation approach. Moreover, the conditions of Routh-Hurwitz [20] have been used to ensure that steady-state solutions are stable and to evaluate their different stability regions.
On the other hand, the employment of the absorbers in the configuration of various constructions of dynamical models has piqued the interest of many researchers, e.g., [21][22][23][24][25][26] due to their applications in the vibrational control of engineering industries. In [21], the authors investigated how a longitudinal absorber may stabilize and regulate the vibration of a ship roll motion through the examination of the behavior of 3-DOF nonlinear spring pendulum. The behavior of a 2-DOF tuned absorber has been studied in [22] when the rotation of its suspension point is on an elliptic trajectory. The frequency response's equations are used to analyze the steady-state solutions around the chosen resonance situation, in which the conditions of stability have been established. The influence of the existence of a damped nonlinear EP on the behavior of this problem is regulated in [23]. The stability analysis for the motion of a transverse absorber linked with a nonlinear elastic spring, according to the examined three cases of the system's resonance, is presented in [24]. The vibrational analysis problem of a connected inverted pendulum with a passive mass for a spring absorber is addressed in [25]. In [26], the authors demonstrated how to modify automatically the rotating speed of a 2-DOF pendulum absorber by determining the phase between the primary structure's vibration and the absorber vibration, in which the response of absorber's speed has been implemented.
Furthermore, there are other types of pendulums that have been studied in [27][28][29][30][31][32]. The motion of a double pendulum under the influence of a vibrating force and a gravity one in the situation of a vibrating point of suspension is investigated in [27]. The parametrically driven double pendulum and its bifurcation structure at tiny oscillation amplitudes are investigated in [28], while the resonance and nonresonance cases of a double pendulum are discussed in [29]. However, the motion of a system comprising of two attached physical pendulums swinging along horizontal axes is investigated in [30]. Examinations of the experimental and numerical results for the motion of a triple pendulum were investigated in [31][32][33][34]. In [31], numerical simulation of generic mechanical model with rigid motion limiters as well as the stability and bifurcation analysis is presented. On the other hand, the experimental rig and its corresponding dynamical model of the triple pendulum are provided in [32]. The piston-connecting rod-crankshaft system of a single-cylinder combustion engine is described in [33] as a particular case of a triple pendulum with barriers. High consistency between the experimental results and the numerical one for the motion of the same model is presented in [34]. The bifurcations and stability of high-frequency periodic motions with limited amplitude are also addressed, as well as the stability criteria. Recently, the stability analysis of un-stretched double pendulum and triple one when they are subjected to external harmonic moments and force are investigated in [35] and [36], respectively. The stability and instability zones for different parameters of the frequency responses are estimated.
This paper focuses on the motion of a plane dynamical system consisting of a TRBP with 3 DOFs under the impact of three harmonic external moments as a novel dynamical model. The governing system of motion, which is made up of a set of three secondorder nonlinear differential equations, is derived using Lagrange's equations. The analytic approximation solutions strategy of this system, as novel solutions, is based on the use of the AMS. The prerequisites for solvability are met by removing secular terms. As a result, the arising resonance cases are identified, and the equations modulation are obtained. The frequency responses curves and the time histories of the derived solutions are drawn to highlight the good effects of various parameters on the motion. The fourth-order Runge-Kutta method is utilized to obtain the numerical results of the governing system. The juxtaposition of the analytical and numerical results demonstrates that the perturbation technique used is extremely accurate. The stability and instability domains are evaluated using Routh-Hurwitz conditions, and the steady-state solutions are analyzed.

Dynamical modeling
The investigated dynamical model in this work consists of three rigid bodies attached with each other and moving in a vertical plane under the influence of a gravity field with acceleration g as shown in Fig. 1.
The physical quantities e j ðj ¼ 1; 2; 3Þ and l k ðk ¼ 1; 2Þ refer to the distance between the bodies gravity centers z j and the rotation centers O j , and the link lengths of the first two pendulums, respectively. Furthermore, the mass and the moments of inertia of jth links regarding to z j (orthogonal to the movement xy plane) are represented by m j and J j , respectively. The excitation moments operating on the links of the pendulums at O j are denoted by M j ðtÞ ¼ f j cos X j t (where f j and X j are the amplitudes and frequencies of these moments), and the angular position's variables are denoted by w j ðj ¼ 1; 2; 3Þ. It is supposed that the resistance moments at O j have a viscous damping with coefficients c j .
Based on the above description, the kinetic energy T and the potential one V of the system can be stated as follows The foregoing set of Eq. (5) provides three secondorder nonlinear differential equations in terms of w j ðj ¼ 1; 2; 3Þ.

The recommended method
The goal of the present section is to use the AMS to acquire the AS of the preceding governing EOM and to classify the different resonances cases. To achieve this aim, the trigonometric function of w j ðj ¼ 1; 2; 3Þ can be approximated in the neighborhood of static equilibrium position up to the third order. Therefore, the system of EOM (5) can be rewritten as follows It is crucial to realize that the vibration amplitudes; specially the generalized coordinates must be formulated in terms of a tiny parameter 0\ e\\1. Therefore, new variables h; /; and v can be introduced as follows w 1 ¼ ehðs; eÞ; w 2 ¼ e/ðs; eÞ; Now, we can seek for these variables as power series of e according to the following series [6] hðs; eÞ ¼ X 3 k¼1 e kÀ1 h k ðs 0 ; s 1 ; s 2 Þ þ Oðe 3 Þ; Taking into account the smallness of the parameters B j ; N j ; M j ; S 1 ; S 2 ; f j ; and C j as below where the parametersB j ;Ñ j ;M j ;S 1 ;S 2 ;f j ;l 1 ;l 2 ; and C j are unity-order parameters. Making use of (11)-(14) into (8)- (10), and then equaling the coefficients of equal powers of e to construct the sets of the following partial differential equations (PDE): Order of e: Order of e 2 : Order of e 3 : The preceding Eqs. (15)-(23) constitute a system of nine linear PDE, which can be solved in the order they appeared. Therefore, the first group of Eqs. (15)- (17) have the general solutions that can be represented in the following ways Here A j ðj ¼ 1; 2; 3Þ and A j are unknown complex functions and the corresponding complex conjugate, respectively.
Inserting (24)-(26) into (18)- (20) and canceling the terms that lead to secular ones to get the following conditions As a result, the solutions of second order can be written as follows where cc denotes the preceding terms' complex conjugates.
Making use of the solutions (24)- (26) and (30)-(32) into Eqs. (21)-(23), one obtains the requirements for deleting secular terms from the approximation of order three in the forms Thus, the third-order approximations h 3 ; / 3 ; and v 3 have the forms The unknown functions A j ðj ¼ 1; 2; 3Þ can be determined by looking at the conditions (27)-(29) and (33)- (35). Therefore, the required AS can be obtained easily after substituting the solutions (24)-(26), (30)-(32), (36)- (38), and series (12) into hypothesis (11). As previously stated, dealing with the AMS necessitates the employment of an unlimited number of separate timescales rather than a single time variable. The solvability constraints, which demand the removal of secular elements from fast time variables, are regarded as the appropriate tax for this flexibility. The fast timescale constrains the structure of the estimated solutions. It is also crucial to double-check that the requirements of the solvability for different orders are consistent. The amplitudes of ''free'' resonant phrases that emerge in each expansion order are limited by these requirements. The analysis may provide inaccurate results or only allow simple solutions if the restrictions are not fulfilled. It should be noted that alternative free amplitude options may yield contradicting findings [38].

Resonance requirements and equations of modulation
The categorizations of resonance cases that may emerge in the second or third-order solutions, as well as the consideration three of these cases, are both significant aspects of this section. It is well known that the resonance cases occur if the denominators of these parts tend to zero [39]. Therefore, it might be categorized as: (i) When p 1 % 1; p 2 % x; and p 3 % -, there is a main primary external resonance. (ii) When x % 1; -% x; and -% 1, there is an internal resonance.
It is noted that, if any one of the resonance cases is met, we should expect that the studied model's behavior will be challenging. Furthermore, the approach described above is valid if the vibrations have values other than that of the resonance. Therefore, we will investigate the three primary external resonances that occur at the same instant to remedy this problem. To pursuit of this objective, it is critical to employ the dimensionless detuning parameters r j ðj ¼ 1; 2; 3Þ that detect the distance between the oscillations and the stern resonance. Then we will be able to write Then, we can formulate these parameters in terms of e according to Inserting (39) and (40) into (18)- (23) and deleting the terms that yield secular ones, to gain the conditions of solvability as follows.
-Conditions of the third-order approximation 2i It is clear that the functions A j ðj ¼ 1; 2; 3Þ can be determined using the forgoing criteria (41)-(46) which can be represented in the following polar form [40] A j ðs 1 ; s 2 Þ ¼ 1 2ã j ðs 1 ; s 2 Þe iw j ðs 1 ;s 2 Þ ; a j ¼ eã j ; h j ðs 1 ; s 2 Þ ¼r j s 1 À w j ðs 1 ; s 2 Þ; where h j are the modified phases, w j are the phase angles and a j are the amplitudes. Making use of (47) into (41)- (46), and then distinguishing the real parts and the imaginary ones of the resulted equations to acquire the following modulation equations of six ordinary differential equations from first order where Fig. 7 Consistency between the AS and the NS of the variable w 1 ; w 2 ; and w 3 when x ¼ 0:4690 Referring to the preceding system (48), the time histories of its solutions a j ðj ¼ 1; 2; 3Þ and h j are graphed in Figs. 2 and 3, respectively. These figures are calculated at xð¼ 0:3829; 0:4062; 0:4690Þ and according to the following data  Fig. 8 to graph the phase planes figures when x takes different values. It is notable that we have closed trajectories which confirm that the gained AS have a stable behavior, which is predicted before, during the tested interval of time.

Steady-state solutions
The goal of this section is to investigate the dynamical system's vibrations under examination the steadystate case. To accomplish this purpose, we consider the null value of the first derivatives of the modified phases h j and amplitudes a j in Eq. (48), i.e., dh j ds ¼ da j ds ¼ 0 ðj ¼ 1; 2; 3Þ [41]. Therefore, the next set of six algebraic equations regarding the variables h j and a j is yielded Fig. 9 Amplitudes' resonance curves a j ðj ¼ 1; 2; 3Þ as a function of r 3 at x ¼ 0:4062 where Elimination of the modified phases h j from the previous system produces three nonlinear algebraic equations of amplitudes a j and the frequency response functions that are clarified by the detuning parameters r j as follows It is important to note that the case of steady-state solutions is considered an important part for the stability's examination. Therefore, to study the behavior near a neighborhood region of fixed points, let us consider the following substitutions into the above system of Eqs. (48) [42,43] a 1 ¼ a 10 þ a 11 ; a 2 ¼ a 20 þ a 21 ; where a j0 ; h j0 are the solutions at the steady state of (49) and a j1 ; h j1 are the corresponding tiny perturbations. Therefore, the linearized system of (48) has the form The solutions of the previous system can be achieved if we represent a j1 and h j1 exponentially as in the forms q k e kT ; where q k ðk ¼ 1; 2; :::; 6Þ are constants and k denotes the unknown perturbations' associated with their eigenvalues. The roots' real parts of the following characteristic equations should not be positive value, if the solutions at steady state a j0 and h j0 are stable asymptotically [44] Here C k are functions of the unperturbed parameters a j0 ; h j0 ; and g j (see Appendix1).
The Routh-Hurwitz criteria [20] provide the necessary and sufficient requirements for steady-state solutions that can be expressed as follows 6 Analysis of system's stability In the present part of this section, the dynamical motion of the considered TRBP system in the presence of external moments M j ðj ¼ 1; 2; 3Þ is examined using the linear stability analysis. The stability requirements are performed besides the modeling of the nonlinear system's equations. Some parameters, such as detuning parameters r j and the natural frequencies 1; -, have been discovered to play a key role in undermining the stability requirements. A customized process with different system parameters has been used to plot the stability graphs of the system (48). The amplitudes a j of the fluctuations versus time are plotted for distinct parametric regions, see Figs. 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, and 23 that are drawn to show the impact of the r j values on the possible fixed points.
Overall, the areas of stable and unstable fixed points under the present constraints of the selected values of the impacted parameters, lie in the range r j À 0:07 and À0:07\r j , respectively. It is important to note that the solid curves reflect the domain of stable fixed points, whereas the dashed curves denote the unstable ranges. No peaks are observed in the planes a 1 r 3 and a 2 r 3 except for three critical fixed points, as shown in parts (a) and (b) of Fig. 9, while there exist one peak and one critical fixed points in the plane a 3 r 3 when x ¼ 0:4062 at ðr 1 ¼ À0:1; r 3 ¼ 0:2Þ; ðr 1 ¼ r 3 ¼ 0Þ; and ðr 1 ¼ 0:1; r 3 ¼ À0:2Þ. It is clear from the third equation in the system (50) that a 3 is a function of r 3 , so we get a curve as shown in Fig. 9c, this means that a 3 changes with r 3 , while we get straight lines as shown in Fig. 9a and b. This means that the value of a 1 , a 2 are invariant with respect to r 3 because they are functions of r 1 and r 2 as observed in the first and second equations of the system (50), respectively. A closer look at Fig. 10 reveals that it is graphed when ðr 1 ¼ 0:1; r 3 ¼ À0:002Þ; ðr 1 ¼ À0:1; r 3 ¼ 0:002Þ; and ðr 1 ¼ r 3 ¼ 0Þ. Several critical points are drawn in part (a). The drawn curve in part (b) includes one peak at the point ð0:1287; 0:0253Þ and one critical fixed point in the plane a 2 r 2 at x ¼ 0:4062: On the other hand, there exist three critical fixed points in Fig. 10c at the same considered values of the detuning parameters.
Curves of Fig. 11 are plotted when ðr 2 ¼ 0:02; r 3 ¼ À0:003Þ; ðr 2 ¼ À0:02; r 3 ¼ 0:003Þ; and ðr 2 ¼ r 3 ¼ 0Þ, the peak fixed point appears in the plane a 1 r 1 of frequency response as in part (a). There are three critical fixed points are observed in the planes of frequency responses a 2 r 1 and a 3 r 1 as shown in Fig. 11b and c.
No peaks are observed in the planes a 1 r 3 and a 2 r 3 except for two and one critical fixed points, as shown in portions (a) and (b), respectively, of Fig. 12. In Fig. 13, there is one critical fixed point in the plane a 1 r 3 and a 2 r 3 while there exist one peak and one critical fixed points in the plane a 3 r 3 of Figs There are two critical fixed points in the plane a 1 r 3 of Fig. 18 while there exists one critical fixed point in part (b) and one peak and one critical fixed points in part (c) at x ¼ 0:4690; C 1 ð¼ 0:005; 0:01; 0:1Þ. In Fig. 19, the one peak and one critical fixed points exist in the plane a 3 r 3 while there is one critical fixed point in the plane a 2 r 3 and a 1 r 3 when x ¼ 0:4690; C 2 ð¼ À0:002; 0:004; 0:006Þ: There are one and two critical points in parts (a) and (b) of Fig. 20, respectively. There exists two critical points and one peak in part (c) at x ¼ 0:4690 and C 3 ð¼ 0:004; 0:006; 0:008Þ.
No peaks are observed in the planes a 1 r 3 and a 2 r 3 except for three critical fixed points, as drawn in parts (a) and (b) of Fig. 21, while there exists essential one peak and one critical fixed points in the plane a 3 r 3 when x ¼ 0:4690 at ðr 1 ¼ À0:1; r 3 ¼ 0:2Þ; ðr 1 ¼ r 3 ¼ 0Þ; and ðr 1 ¼ 0:1; r 3 ¼ À0:2Þ. A closer examination of the curves of Fig. 22 shows that it is graphed when ðr 1 ¼ 0:1; r 3 ¼ À0:002Þ; ðr 1 ¼ À0:1; r 3 ¼ 0:002Þ; and ðr 1 ¼ r 3 ¼ 0Þ. Several critical are drawn in part (a). The drawn curve in part (b) includes two peaks and one critical fixed point in the plane a 2 r 2 at x ¼ 0:4690. On the other hand, there are three crucial critical fixed points in Fig. 22c at the same considered values of the detuning parameters.
Curves of Fig. 23 are drawn when ðr 2 ¼ 0:02; r 3 ¼ À0:003Þ; ðr 2 ¼ À0:02; r 3 ¼ 0:003Þ; and ðr 2 ¼ r 3 ¼ 0Þ, the peak fixed point appears in the plane a 1 r 1 of frequency response as in part (a). There are three observed fundamental critical fixed points in the planes of frequency responses a 2 r 1 and a 3 r 1 as shown in Fig. 23b and c. Tables 1 and 5 show the relation between the critical and peaks fixed points and the values of the parameters r j ; and x, while Tables 2, 3, 4 show the relation between the critical and peaks fixed points at the considered values of the parameters C j and x, respectively (Table 5).
third-order of approximation is investigated using the MSA. The prerequisites for solvability are met by removing secular terms. The emerging resonance cases are detected, and the equations of modulation are achieved. The frequency responses curves and time histories of the obtained new results are visually displayed to emphasize the beneficial effects of the various parameters on the motion. The numerical results of the motion's governing system are gained using the Runge-Kutta fourth-order method. The comparison between both results demonstrates that the used perturbation technique is extremely accurate. The steady-state solutions are investigated, and the zones of stability and instability are evaluated using Routh-Hurwitz conditions. This work is significant because of its applications in engineering vibrational control, which gives the researched system a lot of weight.