EUROPEAN JOURNAL OF PURE AND APPLIED MATHEMATICS 2025, Vol. 18, Issue 3, Article Number 6281 ISSN 1307-5543 – ejpam.com Published by New York Business Global Applying LPA to Analyze the 3D Motion of a Symmetric Rigid Body with a Large Displacement of a Center of Mass T. S. Amer1,∗, A. I. Ismail2, W. S. Amer3, H. F. El-Kafly4 1 Mathematics Department, Faculty of Science, Tanta University, Tanta 31527, Egypt 2 Mechanical engineering Department, College of engineering and architecture, Umm Alqura University, Saudi Arabia 3 Mathematics and Computer Science Department, Faculty of Science, Menoufia University, 32511, Egypt 4Higher Institute of Engineering and Technology, Tanta, Egypt Abstract. This study investigates the 3D motion of a rigid body (RB) rotating around a fixed point, focusing on Lagrange’s case while considering the effects of a gyrostatic moment (GM) and a Newtonian force field (NFF). It is notably that the center of mass is displaced very much from the principal dynamic symmetry axis. Utilizing the fundamental principle of angular momentum, the equations of motion (EOM) are formulated and solved applying the large parameter approach (LPA) to determine approximate solutions (AS) for the irrational frequencies’ case. Euler’s angles, which define the body’s orientation at any moment, are explicitly calculated. Additionally, to assess the impact of applied moments on motion stability, we utilize advanced computational tools to generate graphical representations of the achieved solutions and the related Euler’s angles. This work enhances the understanding of RB dynamics in complex motion scenarios, emphasizing the interaction between external forces and GMs in shaping stability and behavior. The study is poised to greatly influence the aerospace industry by advancing our understanding of rotational motion and the dynamics of celestial bodies, with direct applications in the design and operation of spaceships, spacecraft, and satellites. 2020 Mathematics Subject Classifications: 70E17, 70K45, 34E10, 70F99 Key Words and Phrases: Nonlinear dynamics, Euler’s equations, Gyrostatic moment, Pertur- bation methods, Numerical analysis ∗Corresponding author. DOI: https://doi.org/10.29020/nybg.ejpam.v18i3.6281 Email addresses: tarek.saleh@science.tanta.edu.eg (T. S. Amer), aiismail@uqu.edu.sa (A. I. Ismail), drwaelamer@science.menofia.edu.eg (W. S. Amer), heba elkafly@yahoo.com (H. F. El-Kafly) https://www.ejpam.com 1 Copyright: © 2025 The Author(s). (CC BY-NC 4.0) T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 2 of 35 Table 1: List of abbreviations and famous symbols. RB Rigid body GM Gyrostatic moment NFF Newtonian force field EOM Equations of motion LPA Large parameter approach AS Approximate solutions DEs Differential equations OXY Z Fixed frame Oxyz Rotating frame M Mass of the body λi Gyrostatic moment vector ω = (p, q, r) The body′s angular velocity vector I1, I2, I3 Principal moments of inertia ϵ Large param 1. Introduction For over two centuries, the study of a RB dynamics has been a focal point of research due to its diverse real-world applications. Gyroscopes, for example, play a critical role in stabilizing and guiding the motion of aircraft, submarines, and spacecraft by maintaining balance and determining orientation. The difficulty of the problem is influenced by the initial conditions, particularly the initial angular velocity, center of mass position, and prin- cipal moments of inertia. Moreover, external factors-such as gravitational influences from various centers, whether symmetric or asymmetric, and applied external moments-further complicate the analysis, requiring different mathematical and computational approaches tailored to specific scenarios. These issues are described by the Euler-Poisson’s equations [1], which define the body’s angular velocity and orientation at any given time. Extensive research has explored this topic from various perspectives [2–17]. In [2], the motion of the RB subjected to a combi- nation of gyroscopic and axisymmetric potential forces is analyzed, with two integrals of the system provided. In [3], the authors presented a fourth integral associated with the RB motion around a fixed point under the influence of a gravitational moment. Exact solutions to similar problems necessitate a fundamental fourth first integral, which is only derived in specific well-known cases, such as [2–7]. These cases depend on the location of the body’s center of mass and the values of its main moments of inertia. In [8] and [9], the authors propose specific conditions for solving the equations of Euler-Poisson and simplify them by transforming the system into three nonlinear differential equations (DEs) involving angular velocity components. In [10], an innovative method is introduced for solving the Euler-Poisson system, which is well known for its analytical solution challenges. This new method effectively overcomes these difficulties. In [11], the dynamics of a charged the RB rotating about a fixed point T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 3 of 35 with a rotor attached along a principal axis are examined. The study evaluates sufficient conditions for equilibrium instability using the method of linear approximation. In [12], the spatial rotational motion of a charged symmetric the RB about a stationary point is analyzed, resembling Lagrange’s formulation, where the center of mass is displaced slightly from the body’s dynamic symmetry axis. In [13], the rotational motion of asymmetric the RB influenced by stationary torques in its body-fixed frame, directed on the first component of the gravitational moment, is analyzed. The study presents analytical results and simulations for two scenarios, where a constant torque is applied along the middle and minor axes. The results are interpreted through separatrix surfaces, periodic and non- periodic solutions, equilibrium manifolds, and the extreme values of periodic solutions. The studies in [14] and [15] focus on the motion of a rotating gyro about a fixed point, analyzed using Lagrange-like criteria. By applying the averaging method, a simplified averaged system for the EOM is derived, taking into account a small parameter. This approach yields AS for two scenarios involving different perturbing moments. The author in [16] further explored the rotational behavior of a symmetric gravitational moment near Lagrange’s case, considering the effects of a perturbed torque vector, the third projection of a gravitational moment, and a restoring torque vector. In [17], the PSPM is employed to obtain AS for the motion of the RB in a NFF. Perturbation methods are highly effective in obtaining AS for similar problems; refer to monographs [18–33]For more details, the Krylov-Bogoliubov-Mitropolski technique is extensively utilized in [18] to derive AS for the EOM of a RB’s motion in a gravitational field. However, these solutions include singular points, which are addressed in [19]. In [20–25], the averaging approach is utilized to derive an averaged system for the motion of a symmetric RB from different perspectives. This method is particularly effective in solving the system when the body rotates in a uniform field [20, 21], experiences the GM [22, 23], or involves a charged body [24, 25]. Additionally, in [26], this approach is applied when the center of mass is offset slightly from the axis of dynamic symmetry. The LPA plays a crucial role in studying symmetric rigid bodies experiencing signif- icant rotational or translational displacements. By defining a dominant parameter that characterizes the motion’s scale, AS can be developed to describe the system’s core dynam- ics. This method has been widely applied in fields such as celestial mechanics, robotics, and structural mechanics, where large displacements significantly impact system behavior. In [27, 28], the problem of large RB displacements is examined as a fundamental challenge in classical mechanics and engineering, particularly in relation to nonlinear dynamics and stability. The asymptotic approach based on LPA serves as an effective analytical tool when traditional perturbation methods are insufficient due to strong nonlinearities. The presented works in [29–31] center broadly on the field of dynamics, with a par- ticular emphasis on the motion of rigid bodies. It covers the core principles of analytical dynamics, spacecraft attitude dynamics, robotics mechanics, perturbation techniques for finding AS, and numerical methods for solving the governing DEs. In [32, 33], the authors solved Euler’s equations using three novel approaches to analyze the rotational motion of a charged RB about a fixed point under the influence of GM. Whereas in [34], the authors derived the angular velocity solutions for a RB model, similar to the Lagrange’s scenario T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 4 of 35 with constant torque and GM, by decoupling the EOM and expressing them in complex form. In [35], the authors demonstrate the influence of GMs and the electromagnetic field on the body’s motion. In [36], the averaging technique is used to obtain the averaging system of the EOM for a symmetric RB. Notably, this system has been successfully solved for various scenarios, especially when the body moves in a uniform field. The vibrational motions of RBs are studied in several works, e.g. [37–40]. The planar dynamics of an elastic spring connected to a RB are analyzed in [37], leading to a dynamical model with three-degree-of-freedom that follows the same path as its point of suspension. In [38], a particular case is explored where the suspension point remains fixed. The study in [39] examines linear damped spring that carry a RB, with the pivot point tracing an elliptical path. These works employ the method of multiple-scales to obtain the AS of the EOM. Consequently, solvability requirements are derived by eliminating secular terms, and po- tential resonance cases are defined. In [40], the vibrational behavior of a rigid body around its equilibrium position is investigated numerically using the ode45 solver, based on the fourth-order Runge-Kutta method [41]. The comparison between numerical solutions and asymptotic solutions for the pendulum motions of RBs is presented in [42], demonstrating their reliability and the high accuracy of the perturbation approach utilized. A general- ization of this problem is discussed in [43]. Additionally, the Routh-Hurwitz criteria [44] are applied to confirm the stability of steady-state solutions and to assess their various stability regions. The presented works in [45–48] demonstrate innovative approaches to solving intricate mechanical systems, enhancing the precision and applicability of mathematical models in physics and engineering. In [45], a comparative study of analytical and numerical solutions for a mass-spring system on a massless cart, revealing dynamic interactions that enhance understanding in mechanical engineering contexts is presented. In [46], the authors explored a pendulum situated between two springs using the Ms-DTM method, shedding light on complex dynamic behaviors and the efficacy of this analytical technique. A theoretical study of harmonic oscillators applying the optimal and modified homotopy method is examined in [47] to address time dependence, contributing valuable insights into dynamic system analysis. Moreover, examine a pendulum’s response to vibratory influences using semi-analytical solutions formed through the aforementioned method is examined in [48]. A comprehension of pendulum dynamics and validating sophisticated perturbation techniques within mechanical problems is presented. The work explores the 3D motion of a symmetric the RB with a very much displace- ment of its center of mass from the axis of dynamic symmetry. The system is subjected to the GM and a NFF. The LPA is employed to obtain AS for the EOM in cases of irrational frequencies. Such solutions are then visualized to highlight the influence of various body parameters on the motion over time, serving as a generalization of earlier work [16]. The use of the LPA enables the derivation of AS remain valid even for significant deviations from equilibrium. Recent advancements in asymptotic methods have proven effective in addressing highly nonlinear the RB problems. Euler’s angles for the motion are deter- mined, and graphical representations provide insight into the motion at any given instant. The significance of this work explores the mathematical formulation, derivation, and ap- T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 5 of 35 plications of the LPA in solving large displacement problems for symmetric the RBs. By leveraging asymptotic approximations, we aim to provide practical and computationally efficient solutions to highly nonlinear the RB problems, contributing to both theoretical and applied mechanics. The RB dynamics play a significant role in the theory of gyroscopes, particularly for symmetric RBs that rotate rapidly around their axis of symmetry. Gyroscopes are crucial in numerous engineering applications. For instance, in racing vehicles, the engine itself acts like a large gyroscope, influencing performance and stability. In aviation, gyroscopes help maintain stability and control orientation; when an aircraft turns right or left, the gyroscopic effect can cause the nose or tail to rise or dip accordingly. Similarly, gyroscopes are essential for steering and pitching in ships. There are various types of gyroscopes used to measure angular velocity, including optical, vibrating, and mechanical gyroscopes. Mi- cromechanical gyroscopes are widely used in consumer electronics, as well as in compasses, computer pointing devices, and other systems requiring orientation detection. 2. Overview of the problem This section aims to offer a deeper insight into the dynamic behavior of a heavy symmetric RB with massM , such as the conditions of Lagrange’s gyroscope are satisfied. The body is anchored at a point O, which serves as the origin for two coordinate systems: a stationary one OXY Z and a rotating one Oxyz that is rigidly attached to the body, with its axes aligned along the inertia’s main axes at O. We will suppose that the body is influenced by a GM λ, whose components (λj ; j = 1, 2, 3) are directed on the principal axes Ox, Oy and OZ. Additionally,in which λ1 = 0 and a NFF originating from an attracting center O1, located along the downward fixed axis OZ and far away the origin O, at a significant distance R (see Fig. (1)). Figure 1: The RB’s problem. T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 6 of 35 The EOM are formulated according to the equations Euler-Poisson’s [1], as follows I1ṗ− (I2 − I3)qr + qλ3 − rλ2 =Mg(γyG − βzG)−N(I2 − I3)βγ, I2q̇ − (I3 − I1)rp− pλ3 =Mg(αzG − γxG)−N(I3 − I1)γα, I3ṙ − (I1 − I2)pq + pλ2 =Mg(βxG − αyG)−N(I1 − I2)αβ, α̇ = rβ − qγ, β̇ = pγ − rα, γ̇ = qα− pβ. (1) It is important to note that the inertia principal moments are denoted by Ij along the principal axes. The over dots indicate derivatives regarding time t, g denotes gravitational acceleration, N = 3g/R and ϵ represent a large parameter. In this context, (p, q, r) , (α, β, γ), and xG, yG, zG) refer to the angular velocity’s components, the unit vector along the vertical Z-axis, and the center of mass coordinates in the principal axis frame, respectively. The given EOM have three first integrals associated with energy, area, and geometric one [1] I1p 2 + I2q 2 + I3r 2 − 2Mg(αxG + βyG + γzG) +N(I1α 2 + I2β 2 + I3γ 2) = I1p 2 0 + I2q 2 0 + I3r 2 0 − 2Mg(α0xG + β0yG + γ0zG) +N(I1α 2 0 + I2β 2 0 + I3γ 2 0), I1pα+ (I2q + λ2)β + (I3r + λ3)γ = I1p0α0 + (I2q0 + λ2)β0 + (I3r0 + λ3)γ0, α2 + β2 + γ2 = 1, (2) where p0, q0, r0, α0, β0 and γ0 represent the values of p, q, r, α, β and γ and t = 0, respec- tively. Taking the following parameters and variables into consideration: √ εp = p1 n , √ εq = q1 n , r = r1 n , √ εα = α1, √ ε β = β1, γ = γ1, √ ε λ2 = I2 n λ́2, λ3 = I3 n λ́3, t = nτ, n = √ I1/MgzG, a = (I1 − I3)/I1, b = I3/I1, I1 = I2 ̸= I3, √ εxG = zG, yG = 0, k = Nan2. (3) Upon substituting (3) into (1) and (2), the following equivalent system is obtained dp1 dτ = aq1r1 + r1λ́2 − q1λ́3 − β1 − kβ1γ1, dq1 dτ = −ap1r1 + p1λ́3 + α1 − γ1 + kα1γ1, εdr1dτ = −(p1λ́2 − 1 bβ1), dα1 dτ = r1β1 − q1γ1, dβ1 dτ = p1γ1 − r1α1, ε dγ1 dτ = (q1α1 − p1β1), (4) p21 + q21 − 2α1 + ε(br21 − 2γ1 − kγ21) = p210 + q210 − 2α10 + ε(br210 − 2γ10 − kγ210), p1α1 + (q1 + λ́2)β1 + ε(br1 + λ́3)γ1 = p10α10 + (q10 + λ́2)β10 + ε(br10 + λ́3)γ10, α2 1 + β21 + εγ21 = ε. (5) 3. The used approach This section aims to apply the LPA [27, 28] to achieve the solutions for the RB under investigation. To achieve this aim, we can use the above first integrals to represent r1andγ1, as follows εr1 = εr10 − u2, εγ1 = ε− u1, (6) T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 7 of 35 where u1 = 1 8(4U1 + ε−1U2 1 ) + ..., u2 = 1 8br10 {r210[4(U2 − U20) + ε−1(U2 1 − U2 10)] + ε−2 br210 (U2 − U20) 2 + ...}, U1 = α2 1 + β21 , U2 = p21 + q21 − 2α1 + U1(k + 1). (7) Here, Us0; (s = 1, 2) stands for the values of Us at t = 0 , whereas the dots represent small terms of higher-order powers of ε−1 . To refine system (4) into a more practical form, we remove r1 and γ1 from (4). By incorporating (6) into (4), yields the following nonlinear system of first-order differential equations εdp2dτ = ελ́2r10(1 + χδ) + ερ1q2 +N1, εdq2dτ = −ερ1p2 +N2, εdα2 dτ = εδλ́2r10 + ερ2β2 +N3, εdβ2 dτ = −ερ2α2 +N4, (8) where p1 = p2 + χα2 + d1, q1 = q2 + χβ2, α1 = (1 + χδ)α2 + δp2 + d2, β1 = (1 + χδ)β2 + δq2, χ = −(1 + k)/± √ (br10 + λ́3)2 + 4(1 + k), δ = [−(br10 + λ́3)± √ (br10 + λ́3)2 + 4(1 + k)]/2(1 + k), d1 = r10/(1 + k + r10λ́3 − ar210), d2 = d1/r10, ρ1 = χ(1− r10δ) + (ar10 − δ(1 + k)− λ́3)(1 + χδ), ρ2 = χ[δ(λ́3 − ar10)− 1] + (r10 + δ(1 + k))(1 + χδ), N1 = χ(u2β1 − q1u1)− (1 + χδ)[(aq1 + λ́2)u2 − ku1β1], N2 = χ(u1p1 − u2α1) + (1 + χδ)[ap1u2 + u1 − kα1u1], N3 = δ[(aq1 + λ́2)u2 − ku1β1]− (u2β1 − u1q1), N4 = −δ[ap1u2 + u1 − ku1β1]− (u1p1 − u2α1). (9) To better understand the distribution of frequencies ρs; (s = 1, 2),we choose a suitable value for r0. Consequently, ρ1/ρ2 may be treated as a rational number, leading to a solution of system (8) that is periodic with a period (T0 = 2πns/ρs; s = 1, 2). Next, we reformulate the problem to identify the τ0(ε) solutions of system (8) with a significantly large value of ε. In this context, ε→ ∞ should approach a solution of period T0 of the generating system (8). Therefore, we introduce the following substitution: ετ = (ε+ µ)T, (10) in which µ = µ(ε) that will be estimated latter. Consequently, the problem under investigation is simplified to finding the solutions for the period T0 of the following system εdp2dT = ε(1 + χδ)λ́2r10 + ερ1q2 + V1, εdq2dT = −ερ1p2 + V2, εdα2 dT = εδλ́2r10 + ερ2β2 + V3, εdβ2 dτ = −ερ2α2 + V4, (11) T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 8 of 35 where εV1 = ε[(1 + χδ)λ́2r10 + ρ1q2]µ+ (ε+ µ)N1, εV2 = −ερ1µp2 + (ε+ µ)N2, εV3 = ε(δλ́2r10 + ρ2β2)µ+ (ε+ µ)N3, εV4 = −ερ2α2µ+ (ε+ µ)N4, Vi = V (1) i + ε(−1)V (2) i + ε−2V 3 i + .... (12) We will now proceed to derive the solutions for system (11). Upon closer examina- tion, we observe that this system constitutes a dynamical system comprised of first-order nonlinear DEs involving the variables p2, q2, γ2 and β2. By differentiating these equations once more, we can generate a new system of second-order DEs, which can be expressed as follows ε(d 2p2 dT 2 + ρ21p2) = ρ1V2 + V̇1, ε(d 2q2 dT 2 + ρ21q2) = (−ρ1V1 + V̇2)− ε(1 + χδ)ρ1λ́2r10, ε(d 2α2 dT 2 + ρ22α2) = (ρ2V4 + V̇3), ε(d 2β2 dT 2 + ρ22β2) = (−ρ2V3 + V̇4)− ερ2δλ́2r10. (13) The aforementioned system comprises four distinct nonlinear DEs, each containing a single variable. With this analysis in mind, we can pursue the solutions of this system as follows p2(T, ε) = A1 cos ρ1T +A2 sin ρ1T + ∑ 1, q2(T, ε) = −A1 sin ρ1T +A2 cos ρ1T − [(1 + χδ)λ́2r10]/ρ1 + ∑ 2, α2(T, ε) = A3 cos ρ2T + ∑ 3, β2(T, ε) = −A3 sin ρ2T − δλ́2r10/ρ2 + ∑ 4, ( ∑ d = ∑∞ n=1 ε −nD (n) d (T ); d = 1, 2, 3, 4), (14) with the following initial conditions p2(0, ε) = A1 = A (0) 1 + a1, q2(0, ε) = A2 − [(1 + χδ)λ́2r10]/ρ1 = A (0) 2 + a2 − [(1 + χδ)λ́2r10]/ρ1, α2(0, ε) = A3 = A (0) 3 + a3, β2(0, ε) = −δλ́2r10/ρ2. (15) In this context,A (0) j and aj ; (j = 1, 2.3) represent the unperturbed and perturbed compo- nents Aj , respectively. It is notably that the solutions of system (8), which have T0-periodic solutions in system (13), exhibit a period of εT1 = (ε + µ)T0 . We will assume that µ = µ0 + a4 has solutions of the form (14) that are periodic, seeking ad = ad(ε); d = 1, 2, 3, 4 as a function of the large parameter ε that approaches zero as ε → ∞. In accordance with the LPA procedure, we modify the initial conditions to correspond with the arbitrary constants in the generating system’s solutions. Utilizing equations (9) and (7), we obtain the following T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 9 of 35 results U (0) 1 = d22 + δ2(A2 1 +A2 2) + (1 + χδ)2A2 3 + 2δd2(A1 cos ρ1T +A2 sin ρ1T ) +2d2(1 + χδ)A3 cos ρ2T + 2δ(1 + χδ)A3[A1 cos(ρ1 − ρ2)T +A2 × sin(ρ1 − ρ2)T ] + 2δ2(1 + χδ)λ́2r10( ρ1+ρ2 ρ1ρ2 )(A1 sin ρ1T −A2 cos ρ1T ) +2δ(1 + χδ)2λ́2r10( ρ1+ρ2 ρ1ρ2 )A3 sin ρ2T + δ2(1 + χδ)2(λ́2r10) 2 ×(ρ1+ρ2 ρ1ρ2 )2, U (0) 2 = d21 + (A2 1 +A2 2) + χ2A2 3 − 2d2 + 2(d1 − δ)(A1 cos ρ1T +A2 sin ρ1T ) +2(d1χ− 1− χδ)A3 cos ρ2T + U (0) 1 (1 + k) + 2χA3[A1 cos(ρ1 − ρ2)T +A2 sin(ρ1 − ρ2)T ] + λ́2r10[ χδ(ρ1+ρ2)+ρ2 ρ1ρ2 ][2(A1 sin ρ1T −A2 cos ρ1T ) +2χλ́2r10 ρ1χδ+ρ2(1−χδ) ρ1ρ2 A3 sin ρ2T + (λ́2r10) 2(χδ(ρ1+ρ2)+ρ2 ρ1ρ2 )2]. (16) The expressions for Nd; (d = 1, 2, 3, 4) in (9) can be reformulated in the form N1 = Λ1q2 + Λ2β2 + Λ3, N2 = −(Λ1p2 + Λ2α2) + Λ4, N3 = Λ5q2 + Λ6β2 + Λ7, N4 = −(Λ5p2 + Λ6α2) + Λ8, (17) where Λ1 = −a(1 + χδ)u2 + χ(δu2 − u1) + kβ(1 + χδ)u1, Λ2 = −χa(1 + χδ)u2 − χ[χu1 − (1 + χδ)u2] + kβ(1 + χδ)2u1, Λ3 = −(1 + χδ)λ́2u2, Λ4 = (1 + χδ)[ad1u2 + u1(1− kd2)] + χ(d1u1 − d2u2), Λ5 = δu2(a− 1) + u1(1− kδ2), Λ6 = χ(δau2 + u1)− (1 + χδ)(u2 − ku1), Λ7 = δλ́2u2, Λ8 = −δ[au2d1 + u1(1− kd2)]− (d1u1 − d2u2). (18) Using the aforementioned equations, the unperturbed notations of Λl; (l = 1, 2, ..., 8)are written as: Λ (0) l = [ Γl1 + Γl2 ( A2 1 +A2 2 ) + Γl3A 2 3 ] + Γl4 (A1 cos ρ1T +A2 sin ρ1T ) + Γl5A3 cos ρ2T + Γl6A3 [A1 cos (ρ1 − ρ2)T +A2 sin (ρ1 − ρ2)T ] + Γl7 (A1 sin ρ1T −A2 cos ρ1T ) + Γl8A3 sin ρ1T (19) Where Γlr represents elements of the square matrix∥ Γlr ∥; (l, r = 1, 2, 3, ..., 8). To determine these elements, substitute U (0) 1 , U (0) 2 into equation (7) to obtain u (0) 1 and u (0) 2 . Consequently, Γ1r,Γ2r, ...,Γ8r can be obtained smoothly (see Appendix I). By examining the solutions in equation (14), it becomes clear that these solutions can be achieved if D (n) d (T ) are calculated. By inserting the solutions (14) into system (11) and equaling the coefficients of corresponding powers of ε−1 on both sides, we derive a system T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 10 of 35 that defines the coefficients D (n) d (T ) in (14) as follows dD (n) 1 dT = ρ1D (n) 2 (T ) + V (n) 1 (T ), dD (n) 2 dT = −ρ1D(n) 1 (T ) + V (n) 2 (T ), dD (n) 3 dT = ρ2D (n) 4 (T ) + V (n) 3 (T ), dD (n) 4 dT = −ρ2D(n) 3 (T ) + V (n) 4 (T ), (20) where the conditions D (n) d (0) = 0, are considered and V (n) d (T ) represent known func- tions, if D (e) d is estimated for e < n. At n = 1: Introducing Vd from (12) into (20), matching the coefficients of similar powers of ε−1 on both sides, and differentiating the resulting equations regarding T , we obtain the following equations that determine D (1) d (T ) d2D (1) 1 dT 2 + ρ21D (1) 1 (T ) = Θ11(T ), d2D (1) 2 dT 2 + ρ21D (1) 2 (T ) = Θ21(T ), d2D (1) 3 dT 2 + ρ22D (1) 3 (T ) = Θ31(T ), d2D (1) 4 dT 2 + ρ22D (1) 4 (T ) = Θ41(T ), (21) where, Θ11(T ) = −2ρ1(ρ1µ0 + Λ (0) 1 )p (0) 2 − (ρ1 + ρ2)Λ (0) 2 α (0) 2 + dΛ (0) 1 dT q (0) 2 + dΛ (0) 2 dT β (0) 2 + ρ1Λ (0) 4 + dΛ (0) 3 dT , Θ21(T ) = −2ρ1(ρ1µ0 + Λ (0) 1 )q (0) 2 − (ρ1 + ρ2)Λ (0) 2 β (0) 2 − dΛ (0) 1 dT p (0) 2 −dΛ (0) 2 dT α (0) 2 + dΛ (0) 4 dT − 2ρ1λ́2r10(1 + χδ)µ0 − λ́2r10(1 + χδ)Λ (0) 1 −δλ́2r10Λ(0) 2 − ρ1Λ (0) 3 , Θ31(T ) = −2ρ2(ρ2µ0 + Λ (0) 6 )α (0) 2 − (ρ1 + ρ2)Λ (0) 5 p (0) 2 + dΛ (0) 5 dT q (0) 2 + dΛ (0) 6 dT β (0) 2 + ρ2Λ (0) 8 + dΛ (0) 7 dT , Θ41(T ) = −2ρ2(ρ2µ0 + Λ (0) 6 )β (0) 2 − (ρ1 + ρ2)q (0) 2 − 2ρ2δλ́2r10µ0 − δλ́2r10Λ (0) 6 −λ́2r10(1 + χδ)Λ (0) 5 − ρ2Λ (0) 7 − dΛ (0) 5 dT p (0) 2 − dΛ (0) 6 dT α (0) 2 + dΛ (0) 8 dT . Knowing that Θd1(u), the solutions of the previous system(21) take the form D (1) 1 (T ) = ρ−1 1 ∫ T 0 Θ11(v) sin ρ1(T − v)dv, D (1) 2 (T ) = ρ−1 1 ∫ T 0 Θ21(v) sin ρ1(T − v)dv, D (1) 3 (T ) = ρ−1 2 ∫ T 0 Θ31(v) sin ρ2(T − v)dv, D (1) 4 (T ) = ρ−1 2 ∫ T 0 Θ41(v) sin ρ2(T − v)dv, (22) T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 11 of 35 where, ρ−1 1 Θ11(T ) = G10 +G11 (A1 cos ρ1T +A2 sin ρ1T ) +G12 (A1 cos ρ1T −A2 sin ρ1T ) +G13 cos 2ρ1T +G14 sin 2ρ1T +G15A3 cos ρ2T +G16A3 [A1 cos(ρ1 − ρ2)T +A2 sin(ρ1 − ρ2)T ] +G17A3 [A1 cos(ρ1 − ρ2)T −A2 sin(ρ1 − ρ2)T ] +G18A3 [A1 cos(ρ1 + ρ2)T +A2 sin(ρ1 + ρ2)T ] −G19A3 [A1 cos(ρ1 + ρ2)T −A2 sin(ρ1 + ρ2)T ] +G110 cos(2ρ1 − ρ2)T +G111 sin(2ρ1 − ρ2)T +G112A 2 3 [A1 cos(ρ1 − 2ρ2)T +A2 sin(ρ1 − 2ρ2)T ] +G113A 2 3 [A1 cos(ρ1 − 2ρ2)T −A2 sin(ρ1 − 2ρ2)T ] +G114 cos 2ρ2T +G115 sin 2ρ2T +G116A3 [A1 sin(ρ1 − ρ2)T −A2 cos(ρ1 − ρ2)T ] +G117A3 [A1 sin(ρ1 − ρ2)T +A2 cos(ρ1 − ρ2)T ] +G118A3 [A1 sin(ρ1 + ρ2)T +A2 cos(ρ1 + ρ2)T ] +G119A3 [A1 sin(ρ1 + ρ2)T −A2 cos(ρ1 + ρ2)T ] +G120A3A2 [cos(ρ1 + ρ2)T + cos(ρ1 − ρ2)T ] +G121 (A1 sin ρ1T −A2 cos ρ1T ) +G122A3 sin 2ρ2T +G123 sin(ρ1 − ρ2)T, ρ−1 1 Θ21(T ) = G20 −G21 (A1 sin ρ1T −A2 cos ρ1T ) +G22 (A1 sin ρ1T +A2 cos ρ1T ) +G23 (A1 cos ρ1T +A2 sin ρ1T ) −G24A3 sin ρ2T −G25A3 [A1 sin(ρ1 − ρ2)T −A2 cos(ρ1 − ρ2)T ] +G26A3 [A1 cos(ρ1 − ρ2)T −A2 sin(ρ1 − ρ2)T ] +G27A3 [A1 cos(ρ1 − ρ2)T +A2 sin(ρ1 − ρ2)T ] +G28A3 [A1 sin(ρ1 − ρ2)T +A2 cos(ρ1 − ρ2)T ] −G29A3 [A1 sin(ρ1 + ρ2)T −A2 cos(ρ1 − ρ2)T ] +G210A3 [A1 sin(ρ1 + ρ2)T +A2 cos(ρ1 + ρ2)T ] +G211A3 [A1 cos(ρ1 + ρ2)T −A2 sin(ρ1 + ρ2)T ] +G212A3 [A1 cos(ρ1 + ρ2)T +A2 sin(ρ1 + ρ2)T ] +G213 sin(2ρ1 − ρ2)T +G214A 2 3 [A1 sin(ρ1 − 2ρ2)T +A2 cos(ρ1 − 2ρ2)T ] +G215A 2 3 [A1 sin(ρ1 − 2ρ2)T −A2 cos(ρ1 − 2ρ2)T ] −G216 sin 2ρ2T +G217A3 cos ρ2T +G218A 2 3 cos 2ρ2T, ρ−1 2 Θ31(T ) = G30 +G31A3 sin ρ2T +G32A3 cos ρ2T +G33A 2 3 sin 2ρ2T +G34A 2 3 cos 2ρ2T, T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 12 of 35 ρ−1 2 Θ41(T ) = G40 −G41 (A1 sin ρ1T −A2 cos ρ1T ) +G42 (A1 sin ρ1T +A2 cos ρ1T ) +G43 (A1 cos ρ1T +A2 sin ρ1T ) −G44 sin 2ρ1T +G45 cos 2ρ1T −G46A3 sin ρ2T −G47A3 [A1 sin(ρ1 − ρ2)T −A2 cos(ρ1 − ρ2)T ] +G48A3 [A1 sin(ρ1 − ρ2)T +A2 cos(ρ1 − ρ2)T ] +G49A3 [A1 cos(ρ1 − ρ2)T +A2 sin(ρ1 − ρ2)T ] +G410A3 [A1 cos(ρ1 − ρ2)T −A2 sin(ρ1 − ρ2)T ] −G411A3 [A1 sin(ρ1 + ρ2)T −A2 cos(ρ1 + ρ2)T ] +G412A3 [A1 sin(ρ1 + ρ2)T +A2 cos(ρ1 + ρ2)T ] +G413A3 [A1 cos(ρ1 + ρ2)T +A2 sin(ρ1 + ρ2)T ] +G414A3 [A1 cos(ρ1 + ρ2)T −A2 sin(ρ1 + ρ2)T ] +G415 sin(2ρ1 − ρ2)T +G416 cos(2ρ1 − ρ2)T −G417A 2 3 [A1 sin(ρ1 − 2ρ2)T −A2 cos(ρ1 − 2ρ2)T ] +G418A 2 3 [A1 sin(ρ1 − 2ρ2)T +A2 cos(ρ1 − 2ρ2)T ] −G419 sin(2ρ2)T +G420A3 cos(ρ2T ) +G421 cos(2ρ2)T. Here G1h(h = 0, 1, 2, ..., 23), G2z(z = 0, 1, 2, ..., 18), G3k(k = 0, 1, 2, 3, 4) and (G4u;u = 0, 1, 2, ..., 21) are constants can be acquired easily (see Appendix 2). Substitution of Θd1(u) into (22), yields D (1) 1 (T0) = (E11 +R11)T0, D (1) 2 (T0) = (E21 +R21)T0, D (1) 4 (T0) = −(E31 +R41)T0, (23) where E11 = −1 2 [(G11 −G12)A2 +G121A1], E21 = −1 2 [(G23 −G21)A1 +G22A2], E31 = −1 2G46A3. (24) It is evident that, T0-periodic solutions of (14)must satisfy the following conditions of periodicity for necessity and sufficiency [23]. ψ1 = p2(T0, ε)− p2(0, ε) = 0, ψ2 = q2(T0, ε)− q2(0, ε) = 0, ψ3 = α2(T0, ε)− α2(0, ε) = 0, ψ4 = β2(T0, ε)− β2(0, ε) = 0, (25) Here ψd( d = 1, 2, 3, 4) are treated in terms of (Aj ; j = 1, 2, 3), µ and ε . The conditions in (25) that determine A (0) j , µ0 and aj are not mutually exclusive due to the presence of the first integrals in equations (11) [20]. Furthermore, if A3 ̸= 0 , the third condition is expected to follow directly from the others. It is preferable to treat either A (0) j or α0 as a T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 13 of 35 constant, and express one of the ad as a function of ε which can be removed when ε→ ∞ [31]. Multiplying conditions (25) by ε and setting the free terms of ε to zero, the required conditions for periodicity of D (1) h ( h = 1, 2, 4) are obtained in the form D (1) h (T0) = D (1) h (A1, A2, A3, µ) = 0. (26) To further understand the previous conditions, we utilize equalities (23) to derive: E11 +R11 = 0, E21 +R21 = 0, E31 +R41 = 0. (27) It is evident that if the ratio of the frequencies ρ1 and ρ2 equals 2, 1/2, 1 or −1, the resulting expressions for R11, R21 and R41 are non-zero, as shown below i) for ρ1ρ −1 2 = 2, R11 = −G111 2 , R21 = G216 2 , R41 = A3 [A1 (G48 −G47) +A2 (G49 −G410)] /2. ii) for ρ1ρ −1 2 = 1 2 , R11 = −1 2 [A2A3 (G17 −G16)−A1A3 (G116 +G117)−G123] , R21 = −A3 2 [A1 (G25 −G28) +A2 (G26 −G27)] , R41 = −G44 2 . iii) for ρ1ρ −1 2 = 1, R11 = −1 2 [ A2A 2 3 (G113 −G112) + (G111 +G122A3) ] , R21 = −1 2 [ G213 −G24A3 − (G214 +G215)A1A 2 3 ] , R41 = 1 2 [ (G43 −G41)A1 +G42A2 +G415 + (G417 −G418)A1A 2 3 ] . iv) for ρ1ρ −1 2 = −1, R11 = G122A3 2 , R21 = −G24A3 2 , R41 = 1 2 [A1 (G41 −G43)−G42A2] . (28) Assuming A (0) j ( j = 1, 2, 3) and µ0 satisfy (27), we consider Jacobi’s matrices of Dh(T0)( h = 1, 2, 4) in terms of Aj and compute both of µ for Aj = A (0) j , µ = µ0, and T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 14 of 35 ψd in terms of ad vanishes at ε → ∞ .Computing the second matrix independently of ε . Allowing us to substitute ε → ∞ . Since Aj , µ and ad are carried out in the solutions as related sums, the resulting matrices are identical and can be denoted by J . The solutions of (25) establish a specific case in which periodic solutions exist. 4. The scenario of irrational frequencies This section explores the periodic solutions of the investigated problem when the ratio of ρ1 to ρ2 is irrational. Notably, these solutions arise when A (0) 1 = A (0) 2 = E21 = E31 = 0, E11 ̸= 0, A (0) 3 Γ53 ̸= 0 (where A3 is an arbitrary quantity) [16], J has a third rank, and µ0 = 1 4ρ2 {Γ85 + Γ78 − 2(Γ61 + Γ63A (0)2 3 ) + λ́2r10 ρ2 [δ(Γ65 − 2Γ68) + (1 + χδ)Γ58} (29) are met, and equalities (25) have solutions expressed as power series in terms of the large parameter a1, a2 and a4 where a3 can be set to zero. These solutions approach zero as ε tends to infinity. Based on the above, the periodicity conditions are expressed as follows a1(cos ρ1T0 − 1) + a2 sin ρ1T0 + ε1D (1) 1 (T0) + .... = 0, −a1 sin ρ1T0 + a2(cos ρ1T0 − 1) + ε1D (1) 2 (T0) + .... = 0, D (1) 4 (T0) + .... = 0. (30) According to (29), the unperturbed solutions p20, q20, α20 and β20 have the forms p2(T, 0) = 0, q2(T, 0) = − λ́2r10(1 + χδ) ρ1 , α2(T, 0) = A (0) 3 cos ρ2T, β2(T, 0) = −A(0) 3 sin ρ2T − δλ́2r10 ρ1 . (31) Based on the area integral in (5) and the conditions (15), we obtain A (0) 3 = 1 2G −1 1 [−G2 ± (G2 2 + 4G1G3) 1/2]; G1 ̸= 0, (G2 2 + 4G1G3) > 0. (32) Here Gj(j = 1, 2, 3) are dependently on of h, β, p20, q20, α20 and β20 which can be deter- mined. Utilizing (29) and (30), a1 and a2 are obtained in the forms T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 15 of 35 4εa1 = −T0A(0) 3 {(2 csc2 ρ1 2 T0 − 1)[Γ48 + Γ35 + λ́2r10 ρ1 (δΓ25 − (1 + χδ) ×Γ15)] + 2A (0) 3 [Γ45 + δλ́2r10 ρ2 Γ28 − (1 + 2λ́2r10 ρ1 (1 + χδ)Γ18)]}+ ...., 4εa2 = −T0A(0) 3 {(cot ρ1 2 T0)[Γ48 + Γ35 + λ́2r10 ρ1 (δΓ25 − (1 + χδ)Γ15)] +2A (0) 3 [Γ45 + δλ́2r10 ρ2 Γ28 − (1 + 2λ́2r10 ρ1 (1 + χδ)Γ18)]}+ ..... (33) It is worth noting that, the form of a4 is of order O(ε−2). Utilizing (3), (6), (9), (14), (22), (23), (29), and (33), the desired solutions are expressed as power series of ε−1 in the following forms p = ε−1/2n−1[d1 + χA (0) 3 cos(ρ2n −1t)] + ε−3/2n−1{ χ 2ρ2 [2(Γ81 + Γ83)− Γ65A (0)2 3 ] + 1 2ρ1 [2Γ41 − Γ27 − Γ18 + (2Γ43 − Γ25)A (0)2 3 ] + { χt 2nA (0) 3 [Γ88 − Γ75 + λ́2r10 ρ2 ×(δΓ65 − (1 + δχ)Γ55)]− χ 2ρ2 [2Γ81 + (2Γ83 − Γ65)A (0)2 3 ]− χ 6ρ2 Γ65A (0)2 3 } × cos(ρ2n −1t) + {− χ 2ρ2 A (0)2 3 [Γ88 − Γ75 + λ́2r10 ρ2 (δΓ65 − (1 + χδ)Γ55)] + χt 2nA (0) 3 [−2(Γ61 + Γ63A (0)2 3 )] + Γ85 + Γ75 − 2ρ2µ0 − λ́2r10 ρ2 (δΓ68 − (1 + χδ) ×Γ58) + χ ρ2 Γ68A (0)2 3 ]} sin(ρ2n−1) + {−T0 4 A (0) 3 (2 csc2 ρ1 2 T0 − 1)[Γ48 + Γ35 + λ́2r10 ρ1 (δΓ25 − (1 + χδ)Γ15)]− [ 1 2ρ1 [2Γ41 − Γ27 − Γ18 + (2Γ43 − Γ25)A (0)2 3 ] + t 2nΓ46A (0) 3 ] + T0 2 A (0)2 3 [Γ45 − Γ18 + λ́2r10 ρ1 (δΓ28 − 2(1 + χδ)Γ18)]} × cos(ρ1n −1t) + {1 4 [ 2 ρ1 − T0A (0) 3 (cot ρ1 2 T0)[Γ48 + Γ35 + λ́2r10 ρ1 (δΓ28 −(1 + χδ)Γ15)] + 1 2A (0)2 3 T0[Γ45 + δλ́2r10 ρ2 Γ28 − (1 + 2λ́2r10 ρ1 (1 + χδ)Γ18)] +1 2A (0)2 3 T0[Γ45 − Γ18 + λ́2r10 ρ1 (δΓ28 − 2(1 + χδ)Γ18)]− t 6nρ1 A (0)2 3 T0Γ28 +A (0) 3 T0[Γ45 − Γ38 − λ́2r10 ρ1 (δΓ28 − (1 + χδ))]} sin(ρ1n−1t) + χ 6ρ2 [A (0)2 3 Γ65 × cos(2ρ2n −1t)− 3Γ68 sin(2ρ2n −1t)] + χ 6ρ1 A (0)2 3 Γ28 sin(2ρ1n −1t)}+ ...., T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 16 of 35 q = −ε−1/2n−1{χ[ δλ́2r10 ρ2 +A (0) 3 sin(ρ2n −1t)]− λ́2r10 ρ1 (1 + χδ)}+ ε−3/2n−1{− 1 ρ2 {χ[2δλ́2r10µ0 − Γ71 + λ́2r10 ρ2 [δ(Γ61 + Γ63A (0)2 3 )− (1 + χδ)(Γ51 + Γ53A (0)2 3 )] +1 2A (0)2 3 (Γ68 − 2Γ73)] + 1 ρ1 [−2(1 + χδ)λ́2r10µ0 − Γ31 + 1 2A (0)2 3 (Γ28 − 2Γ33)] − λ́2r10 ρ1 [δ(Γ21 + Γ23A (0)2 3 )− (1 + χδ)(Γ11 + Γ13A (0)2 3 )]}+ χ{− t 2nA (0) 3 [Γ85 + Γ78 −4ρ2µ0 − 2(Γ61 + Γ63A (0)2 3 ) + λ́2r10 ρ2 [δ(Γ65 − 2Γ68) + (1 + χδ)Γ58]] + 1 ρ2 [δλ́2r10µ0 −(Γ71 + Γ73A (0)2 3 ) + λ́2r10 ρ2 [δ(Γ61 + Γ63A (0)2 3 )− (1 + χδ)(Γ51 + Γ53A (0)2 3 )]]} × cos(ρ2n −1t) + χ 2A (0) 3 { t n [Γ88 − Γ75 + λ́2r10 ρ2 [δΓ65 − (1 + χδ)Γ55]]− 2 ρ2 Γ65A (0) 3 + 1 ρ2 [Γ85 + Γ78 − 4ρ2µ0 − 2(Γ61 + Γ63A (0)2 3 ) + λ́2r10 ρ2 [δ(Γ65 − 2Γ68) + (1 + χδ)Γ58]]} × sin(ρ2n −1t) + χ 2ρ2 A (0)2 3 [Γ68 cos(2ρ2n −1t) + Γ65 sin(2ρ2n −1t)] + T0 4 A (0) 3 {(2 × csc2 ρ1 2 T0 − 1)[Γ48 + Γ35 + λ́2r10 ρ1 (δΓ25 − (1 + χδ)Γ15)]− T0 2 A (0) 3 [Γ45 − Γ18 + λ́2r10 ρ1 (δΓ28 − 2(1 + χδ)Γ18)]− 1 2ρ1 [Γ18 − Γ45 + λ́2r10 ρ1 (δΓ28 + 2(1 + χδ)Γ18)] − 1 3ρ1 Γ25A (0) 3 + t 2n [Γ15 + Γ18 + λ́2r10 ρ1 (δΓ25 + 2(1 + χδ)Γ15)]} sin(2ρ1n−1t) +{−T0 4 A (0) 3 cot(ρ12 T0)[Γ48 + Γ35 + λ́2r10 ρ1 (δΓ25 − (1 + χδ)Γ15)] + 1 2A (0)2 3 T0[Γ45 −Γ18 + λ́2r10 ρ1 (δΓ28 − 2(1 + χδ)Γ18)]− 1 ρ2 [−2(1 + χδ)λ́2r10µ0 + 1 2Γ28A (0)2 3 −(Γ31 + Γ33A (0)2 3 )]− λ́2r10 ρ1 [δ(Γ21 + Γ23A (0)2 3 )− (1 + χδ)(Γ11 + Γ13A (0)2 3 )] + t 2nA (0) 3 [Γ18 − Γ45 + λ́2r10 ρ1 (δΓ28 + 2(1 + χδ)Γ18)] + 1 6ρ1 Γ28A (0)2 3 )} cos(ρ1n−1t) + 1 6ρ1 A (0)2 3 [Γ25 sin(2ρ1n −1t)− Γ28 cos(2ρ1n −1t)]}+ ..., r = r0 n − ε−1 nbr10 A (0) 3 {[(1 + χδ)(d2 − 1) + χd1](cos(ρ2n −1t)− 1) + λ́2r10 ρ1ρ2 [δ(ρ1 + ρ2)[χ 2 + (1 + χδ)2] + ρ2χ] sin(ρ2n −1t)}+ ...., α = ε−1/2[d2 + (1 + χδ)A (0) 3 cos(ρ2n −1t)] + ε−3/2{ (1+χδ) 2ρ2 [2Γ81 + (2Γ83 − Γ65) ×A(0)2 3 ] + δ 2ρ1 [2Γ41 − Γ27 − Γ18 + (2Γ43 − Γ25)A (0)2 3 ] + (1 + χδ){ t 2nA (0) 3 ×[Γ88 − Γ75 + λ́2r10 ρ2 (δΓ65 − (1 + χδ)Γ55)]− 1 2ρ2 [2Γ81 + (2Γ83 − Γ65)A (0)2 3 ] − 1 6ρ2 A (0)2 3 Γ65} cos(ρ2n−1t) + (1 + χδ){− 1 2ρ2 A (0) 3 [Γ88 − Γ75 + λ́2r10 ρ2 (δΓ65 −(1 + χδ)Γ55)] + t 2nA (0) 3 [Γ85 + Γ78 − 2ρ2µ0 − 2(Γ61 + Γ63A (0)2 3 − λ́2r10 ρ2 (δΓ68 −(1 + χδ)Γ58)] + 1 ρ2 A (0)2 3 Γ68} sin(ρ2n−1t) + (1+χδ) 6ρ2 A (0)2 3 [Γ65 cos(2ρ2n −1t) −Γ68 sin(ρ2n −1t)]− {T0δ 4 A (0) 3 (2 csc2 ρ1 2 T0 − 1)[Γ48 + Γ35 + λ́2r10 ρ1 (δΓ25 −(1 + χδ)Γ15)]− T0 2 A (0)2 3 [Γ45 − Γ18 + λ́2r10 ρ1 (δΓ28 − 2(1 + χδ)Γ18)] + δ 2 [ 1 ρ1 (2Γ41 −Γ27 − Γ18) + 1 ρ1 (2Γ43 − Γ25)A (0)2 3 + t nΓ46A (0) 3 ]} cos(ρ1n−1t)− δ{−1 4 [ 2 ρ1 −T0A(0) 3 (cot ρ1 2 T0)[Γ48 + Γ35 + λ́2r10 ρ1 (δΓ25 − (1 + χδ)Γ15)] + 1 2A (0)2 3 T0[Γ45 −Γ18 + λ́2r10 ρ1 (δΓ28 − 2(1 + χδ)Γ18)] + t 6ρ1 Γ28A (0)2 3 +A (0) 3 [Γ45 − Γ38 + λ́2r10 ρ1 ×(1 + χδ − δΓ28)]} sin(ρ1n−1t) + δ 6ρ1 A (0)2 3 sin(2ρ1n −1t)}+ ...., T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 17 of 35 β = −ε−1/2(1 + χδ){ δλ́2r10 ρ1ρ2 (ρ1 + ρ2) +A (0) 3 sin(ρ2n −1t)}+ ε−3/2{− (1+χδ) 2ρ2 [2δλ́2r10 ×µ0 + λ́2r10 ρ2 [δ(Γ61 + Γ63A (0)2 3 )− (1 + χδ)(Γ51 + Γ53A (0)2 3 )]− (Γ71 + Γ73 ×A(0)2 3 ) + 1 2A (0)2 3 Γ68] + δ ρ1 [−2(1 + χδ)λ́2r10µ0 + 1 2A (0)2 3 (Γ28 − 2Γ33)− Γ31 − λ́2r10 ρ1 [δ(Γ21 + Γ23A (0)2 3 )− (1 + χδ)(Γ11 + Γ13A (0) 3 )]] + (1 + χδ){− t 2n}[Γ85 +Γ78 − 4ρ2µ0 + λ́2r10 ρ2 [δ(Γ65 − 2Γ68) + (1 + χδ)Γ58]] + 1 ρ2 [2δλ́2r10µ0 − Γ71 −Γ73A (0)2 3 + λ́2r10 ρ2 [δ(Γ61 − Γ63A (0)2 3 )− (1 + χδ)(Γ51 + Γ53A (0)2 3 )]]} × cos(ρ2n −1t) + (1 + χδ){ t 2nA (0) 3 [Γ88 − Γ75 + λ́2r10 ρ2 [δΓ65 − (1 + χδ)Γ55]]− 1 ρ2 ×Γ65A (0)2 3 + 1 2ρ2 A (0) 3 [Γ85 + Γ78 − 4ρ2µ0 + λ́2r10 ρ2 [δ(Γ65 − 2Γ68) + (1 + χδ)Γ58 −2(Γ61 + Γ63A (0)2 3 )]]} sin(ρ2n−1t) + (1+χδ) 2ρ2 A (0)2 3 [Γ68 cos(2ρ2n −1t) + Γ65 × sin(2ρ2n −1t)] + δ{T0 4 A (0) 3 (2 csc2 ρ1 2 T0 − 1)[Γ48 + Γ35 + λ́2r10 ρ1 (δΓ25 − (1 + χδ)Γ15)] +T0 2 A (0)2 3 [Γ18 − Γ45 − λ́2r10 ρ1 (δΓ28 − 2(1 + χδ)Γ18)]− 1 3ρ1 A (0)2 3 Γ25 − 1 2ρ1 A (0) 3 [Γ18 −Γ45 + λ́2r10 ρ1 (δΓ28 + 2(1 + χδ)Γ18)] + t 2nA (0) 3 [Γ18 + Γ15 + λ́2r10 ρ1 (δΓ25 +2(1 + χδ)Γ15)]} sin(ρ1n−1t) + δ{−T0A (0) 3 4 (cot ρ1 2 T0)[Γ48 + Γ35 + λ́2r10 ρ1 (δΓ25 − (1 + χδ)Γ15)] +1 2A (0)2 3 T0[Γ45 − Γ18 + λ́2r10 ρ1 (δΓ28 − 2(1 + χδ)Γ18)] + 1 ρ1 [2λ́2r10(1 + χδ)µ0 − 1 2Γ28 ×A(0)2 3 + Γ31 + Γ33A (0)2 3 + λ́2r10 ρ1 [δ(Γ21 + Γ23A (0)2 3 )− (1 + χδ)(Γ11 + Γ13A (0)2 3 )]] + t 2n ×A(0) 3 [Γ18 − Γ45 + λ́2r10 ρ1 (δΓ28 + 2(1 + χδ)Γ18)] + 1 6ρ1 Γ28A (0)2 3 } cos(ρ1n−1t) + δ 6ρ1 A (0)2 3 [Γ25 sin(2ρ1n −1t)− Γ28 cos(2ρ1n −1t)]}+ ..., γ = 1− 1 2ε −1{(1 + χδ)2[ δλ́2r10 ρ2 +A (0) 3 sin(ρ2n −1t)]2 + [d2 + (1 + χδ)A (0) 3 cos(ρ2n −1t)]2}+ ...., µ0 = 1 4ρ2 {Γ85 + Γ78 − 2(Γ61 + Γ63A (0)2 3 ) + λ́2r10 ρ2 [δ(Γ65 − 2Γ68) + (1 + χδ)Γ58]}, T0 = 2πρ−1 2 . (34) As a result of the fact that the rotatory motion of a RB in 3D is governed by non- linear differential equations, it is frequently difficult or even impossible to obtain exact solutions in many practical scenarios [1, 3, 7]. Therefore, the AS play an essential role in comprehending, predicting, and optimizing the rotational behavior of such systems. The following are the primary reasons why that are important: This problem is governed by Euler’s equations [49], which are highly nonlinear due to the interdependence of angular velocities and moments of inertia. The nonlinear nature of the problem makes finding exact solutions for arbitrary initial conditions or external forces nearly impossible, necessitating the use of approximate methods such as perturbation techniques, numerical simulations, or asymptotic analysis [50]. The use of the AS make it possible to do computations in an effective manner, which makes them suited for real-time applications in which exact results are not feasible. For the T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 18 of 35 purpose of providing feedback control, numerous control systems in the fields of robotics, aircraft, and navigation require extremely quick and real-time estimates of rotational mo- tion [51]. By analyzing small perturbations, AS facilitate the assessment of rotational stabil- ity, enabling engineers to predict whether a system will sustain stable motion or shift into chaotic behavior. The evolution of rotational motion in RBs over extended periods, including effects like precession and nutation, can be computationally demanding [52]. The AS provide an efficient alternative for studying these dynamics. Rotating systems are often susceptible to resonance, where periodic external forces intensify oscillations. Approximate approaches like multiple-scales method and averaging one [50] help detect critical frequencies and optimize system design to prevent instability. RBs in engineering applications frequently undergo compound rotational motion, in- cluding spinning and chaotic dynamics. AS provide computationally efficient models that retain essential physical insights. Numerical simulations, including finite element analysis and Runge-Kutta integration, provide precise results but are computationally demanding [53]. Approximate analytical methods offer initial insights, serving as a foundation for further refinement using computational techniques. As benchmark models, the AS help validate experimental data, allowing researchers to compare predictions, refine models, and detect discrepancies. Therefore, the nonlinear nature of the governing EOM makes obtaining AS RB ro- tation essential. These approximations enable real-time computation, stability analysis, and long-term motion prediction, providing practical benefits across aerospace, robotics, biomechanics, and structural engineering. 5. Euler’s angles and interpretation of motion This section is dedicated to analyzing Euler’s angles for the motion of the given body within the context of the obtained solutions, offering a geometrical interpretation of motion through at any moment. These angles are crucial for determining the orientation of a RB in 3D space. They offer a structured approach to defining a RB’s orientation in 3D space through three sequential rotations relative to a fixed reference frame, ensur- ing a clear and systematic representation of rotation [54]. With only three parameters, these angles provide a simpler and more computationally efficient alternative to rotation matrices and quaternions, making them easier to interpret in engineering contexts. Euler angles are essential in aerospace for defining the orientation of aircraft, satellites, and spacecraft, supporting navigation, attitude control, and stabilization. In robotics, they facilitate precise positioning and motion planning for robotic arms and autonomous systems [55]. They are fundamental to RB dynamics, aiding in the formulation and solution of Euler’s EOM. They offer a simplified representation of angular velocity and acceleration, enhancing mechanical and aerospace studies. Numerous sensors, such gyroscopes are useful for motion tracking and control systems in the real world because they measure orientation in terms of these angles. For smooth and organic object orientation control, Euler angles are frequently utilized in computer graphics T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 19 of 35 and animation. Additionally, they are employed in simulation and game development physics engines. They provide a structured approach to analyzing rotational behaviors, including precession, nutation, and spin, which are essential in planetary motion, satellite dynamics, and gyroscopic systems [55]. A fundamental approach involves applying the expressions for Euler’s angles θ, ψ and ϕ as in [1]. θ = cos−1 γ, ψ̇ = pα+ qβ 1− γ2 , ϕ̇ = r − ψ̇ cos θ, ϕ0 = tan−1 ( α0 β0 ) ; ( ·̇ ≡ d dt ) . (35) Substituting the achieved solutions (34) into expressions (35) to obtain θ, ψ and ϕ as follows θ = θ0 + ε−1(1 + χδ)A (0) 3 csc θ0 [ δλ′ 2r10 ρ2 (1 + χδ) [ sin ( ρ2n −1(t+ χ) ) − sin ( ρ2n −1χ )] +d2 [ cos ( ρ2n −1(t+ χ) ) − cos ( ρ2n −1χ )]] + · · · , ψ = ψ0 + εn−1 csc2 θ0 { d1d2 + (1 + χδ) [ δ(λ́2r10)2 ρ21ρ2 (ρ1 + ρ2) [ δχ ρ2 (ρ1 + ρ2) + 1 ] + χA (0)2 3 ] t + nA (0) 3 ρ2 [d2χ+ d1(1 + χδ)] sin ( ρ2n −1t ) − nλ́2r10 ρ1ρ2 (1 + χδ)A (0) 3 [ 2χδ ρ2 (ρ1 + ρ2) + 1 ] × ( cos ( ρ2n −1t ) − 1 )} + · · · , ϕ = ϕ0 + r10 n t− ε−1 n { εA (0) 3 nbr10 [(1 + χδ)(d2 − 1) + χd1] [ n ρ2 sin ( ρ2n −1t ) − t ] −nλ́2r10 ρ1ρ22 [ δ(ρ1 + ρ2) ( χ2 + (1 + χδ)2 ) + ρ2χ ] ( cos ( ρ2n −1t ) − 1 )} +csc2 θ0 cos θ0 { d21 + (1 + χδ) [ δ(λ́2r10)2 ρ21ρ2 (ρ1 + ρ2) [ χδ ρ2 (ρ1 + ρ2) + 1 ] + χA (0)2 3 ] t + nA (0) 3 ρ2 [d2χ+ d1(1 + χδ)] sin ( ρ2n −1t ) − nλ́2r10 ρ1ρ2 (1 + χδ)A (0) 3 [ 2χδ ρ2 (ρ1 + ρ2) + 1 ] × ( cos ( ρ2n −1t ) − 1 )} + · · · (36) where ϕ0 = π 2 − ρ1n −1χ. Equation (36) clearly shows that the Euler angles of the examined motion are depen- dent on the arbitrary constants θ0, ψ0 and r0. T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 20 of 35 6. Evaluation of the outcomes This section primarily focuses on presenting and analyzing the results outlined in (34) and (36) using the graphical representations of the obtained solutions to highlight the influence of various parameters on the body’s motion over time and to emphasize the significance of the GM. Accordingly, the following data are considered I1 = I2 = 103kg.m2, I3 = 800kg.m2, ε = 1000, xG = 1m, M = 50kg, g = 9.8kg.s−2, r0 = 1.2, θ0 = 10−10 λ2, λ3 = (10, 30, 60)kg.m2.s−1. Fig.2, and Fig.4 illustrate the temporal progression of the outcomes obtained from the system of equations (34) under different values of λ2 and λ3 , respectively. Meanwhile, the corresponding phase plane diagrams, generated using the same parameter settings, are depicted in Figs. Fig.3, and Fig.5. These visualizations are derived from the data previously presented. Fig.2 illustrates the variation of the time-dependent functions p, q, r, α, β and γ when λ2 has distinct values of (10, 30, 60)kg.m 2.s−1 . Certain sections of this figure are computed at λ3 = 10kg.m2.s−1 . The main purpose of the curves in Fig.2 is to examine the positive effects of the selected values of the gyro’s second projection λ2 = (10, 30, 60)kg.m2.s−1 on the body’s behavior. As anticipated, periodic waves appear in certain portions of this figure. The solutions p and α show only small variations with different values of λ2 , as demonstrated in Fig.2 in a and d respectively. This behavior arises from the governing equations, which either do not explicitly depend on it or encompass it within a bracket alongside λ2. However, the other potions (b), (c), (e), and (f) are significantly influenced by changes in values. Notably, the wave amplitudes increase as increases, as shown in Fig.2 in e and f respectively. Additionally, neither the number of fluctuations nor the wavelength of the waves undergoes any change. The body’s steady behavior is depicted through the phase plane graphs of the obtained solutions when λ2(= 10, 30, 60)kg.m2.s−1 , as shown in Fig.3. The plotted curves represent the solutions via their first derivatives, with time being excluded. Notably, symmetric closed curves emerge, highlighting and validating the stability behavior of the body’s motion. Curves in Fig.4 show how the solutions p, q, r, α, β and γ change for distinct values of λ3 (e. g.10, 30 and 60). These solutions are affected by variations in λ3 , producing periodic waves with distinct amplitudes, as shown in parts of Fig.3. As values increase, the amplitudes of these waves decrease. The plots of phase plane of these waves, shown in portions of Fig.5, further demonstrate the stable nature of the results. The results indicate that the solutions maintain stability, without any chaotic characteristics. It is evident that Euler’s angles θ, ψ and ϕ significantly contribute to determining the orientation of the RB. To investigate this aspect, relations (37) are plotted in Figs.6 and Figs.7, considering the previously provided data to analyze the impact of different values of λ2, and λ3 on the body’s dynamic behavior. Notably, Figs.6 are plotted at T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 21 of 35 λ2(= 10, 30, 60)kg.m2.s−1 and λ3 = 10 respectively, while Figs. (7) are computed at λ2 = 10 and,λ3(= 10, 30, 60)kg.m2.s−1 respectively. Notably, the behavior of the nutation angle θ positively influences λ3 as it increases over time, as illustrated in Fig.6 and Fig.7, in panel a . This effect arises from the first approximate term, which depends on ρ1, ρ2, χ and δ . It is apparent that the wave amplitude grows moderately, while the oscillation count remains constant. On the other side, the behavior of the waves representing the precession angle ψ changes gradually as λ2 increases, as shown in panel b of Fig.6 and Fig.7. However, there is no noticeable variation in the waves describing the self-rotation angle ϕ , as seen in panel c of Fig.6 and panel c of Fig.7. The negative values of ϕ arise due to the difference between the first and second terms in the equation, where the second term is significantly larger than the first. Additionally, it should be noted that the inaccuracy in ϕ plays a crucial role in these results. Figure 2: Presents the variations of p(t), q(t), r(t), α(t), β(t) and γ(t) when λ2 = (10, 30, 60) kg ·m2 · s−1 and λ3 = 10 are considered. T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 22 of 35 Figure 3: Depicts the phase plane representations of time dependent functions in Fig.2. Figure 4: Shows the variation of p(t), q(t), r(t), α(t), β(t) and γ(t) when λ3 = (10, 30, 60)kg.m2.s−1 and λ2 = 10 . T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 23 of 35 Figure 5: Reveals the phase plane plots of time dependent functions in Fig.4. Figure 6: Presents the temporal histories of: (a) θ, (b) ψ, and (c) ϕ at λ2(= 10, 30, 60) and λ3 = 10. T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 24 of 35 Figure 7: Shows the time histories of the angles (a) θ, (b) ψ, and (c) ϕ when λ3(= 10, 30, 60) and λ2 = 10 are consided. 7. Conclusion This work analyzes the three-dimensional motion of a new RB according to La- grange’s gyroscope model, considering a new condition by giving the body a very much displacement of the center of mass from the axis of dynamic symmetry. The body is in- fluenced by a GM along the main inertia axes and a NFF. The governing EOM and their first integrals are reformulated as a nonlinear dynamical system of four first-order DEs, which are solved via the LPA. Advanced computational tools are employed to analyze the solutions and visualize the influence of varying applied moments and NFF on the body’s motion. These solutions extend and generalize previous findings by Elfimov [26], which considered the absence of external moments and forces, and by Amer [16], which exam- ined the case of a uniform field. The stability of the obtained solutions is evaluated using computational methods, generating diagrams of the phase plane which confirm the motion remains free of chaotic behavior. Additionally, geometric interpretations are conducted to determine the body’s orientation over time, utilizing Euler’s angles. Time-histories plots of these angles further illustrate the effects of different physical parameters on motion. This research is highly relevant due to its wide-ranging applications in celestial mechanics, physics, and gyroscopic systems, offering valuable insights into the dynamics of rigid bod- ies subjected to complex force interactions. This study is expected to significantly impact the aerospace industry by enhancing our understanding of rotational motion and celestial dynamics, with direct applications in the design and operation of spacecraft, satellites, and other space vehicles. T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 25 of 35 Authors’ Statements • T. S. Amer: Investigation, Methodology, Data curation, Conceptualization, Vali- dation, Reviewing and Editing. • A. I. Ismail: Methodology, Data curation, Conceptualization, Validation, Review- ing and Editing. • W. S. Amer: Conceptualization, Resources, Formal Analysis, Validation, Visual- ization, and Reviewing. • H. F. El-Kafly: Resources, Methodology, Conceptualization, Validation, Formal Analysis, Visualization, and Reviewing. Conflict of Interest The authors declare that they have no conflict of interest. Funding Acknowledgements The author extends their appreciation to Umm-Alqura University, Saudi Arabia for funding this research work through grant number: 25UQU4240002GSSR06 Data Availability Data sharing not applicable to this article as no datasets were generated or analyzed during the current study. References [1] Eugene Leimanis. The General Problem of the Motion of Coupled Rigid Bodies about a Fixed Point, volume 7. Springer Science & Business Media, 2013. [2] H. M. Yehia and A. A. Elmandouh. New conditional integrable cases of motion of a rigid body with kovalevskaya’s configuration. Journal of Physics A: Mathematical and Theoretical, 44(8), 2011. [3] A. I. Ismail and T. S. Amer. A necessary and sufficient condition for solving a rigid body problem. Technische Mechanik, 31(1):50–57, 2011. [4] Iu. A. Arkhangel’skii. On the algebraic integrals in the problem of motion of a rigid body in a newtonian field of force. Journal of Applied Mathematics and Mechanics, 27(1):247–254, 1963. [5] Hamad M. Yehia. Rigid body dynamics. Advances in Mechanics and Mathematics, 2022. [6] A. A. Elmandouh. New integrable problems in a rigid body dynamic with cubic integral in velocities. Results in Physics, 8:559–568, 2018. T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 26 of 35 [7] T. S. Amer and W. S. Amer. The substantial condition for the fourth first integral of the rigid body problem. Mathematics and Mechanics of Solids, 23(8):1237–1246, 2018. [8] S. V. Ershkov. A riccati-type solution of euler-poisson equations of rigid body rotation over the fixed point. Journal of Aerospace Engineering, 228:2719–2723, 2017. [9] S. V. Ershkov and D. Leshchenko. On a new type of solving procedure for euler- poisson equations (rigid body rotation over the fixed point). Acta Mechanica, 230:871– 883, 2019. [10] S. V. Ershkov and D. D. Leshchenko. On a new type of solving procedure for euler-poisson equations (rigid body rotation over the fixed point). Acta Mechanica, 230(1):871–883, 2019. [11] A. A. Elmandouh. On the stability of the permanent rotations of a charged rigid body-gyrostat. Acta Mechanica, 228:3947–3959, 2017. [12] T. S. Amer, I. M. Abady, H. A. Abdo, and H. F. El-Kafly. The asymptotic solutions for the motion of a charged symmetric gyrostat in the irrational frequency case. Scientific Reports, 14:16662, 2024. [13] T. S. Amer, H. F. El-Kafly, A. H. Elneklawy, and A. A. Galal. Analyzing the spa- tial motion of a rigid body subjected to constant body-fixed torques and gyrostatic moment. Scientific Reports, 14(1):5390, 2024. [14] T. S. Amer and I. M. Abady. On the motion of a gyro in the presence of a newtonian force field and applied moments. Mathematics and Mechanics of Solids, 23(9):1263– 1273, 2018. [15] T. S. Amer. On the rotational motion of a gyrostat with mass distribution. Nonlinear Dynamics, 54:189–198, 2008. [16] T. S. Amer. On the motion of a gyrostat similar to lagrange’s gyroscope under the influence of a gyrostatic moment vector. Nonlinear Dynamics, 54:249–262, 2008. [17] M. A. Farag. Analysis of the rotational motion of a solid body in the presence of external moments. Journal of Vibration Engineering & Technologies, 12(1):757–771, 2023. [18] A. I. Ismail. On the application of krylov-bogoliubov-mitropolski technique for treat- ing the motion about a fixed point of a fast spinning heavy solid. Zeitschrift für Flugwissenschaften und Weltraumforschung (ZFW), 20(4):205–208, 1996. [19] A. I. Ismail. Treating a singular case for a motion of rigid body in a newtonian field of force. Archive of Mechanics, 49(6):1091–1101, 1997. [20] D. D. Leshchenko and A. S. Shamaev. Perturbed rotational motions of a rigid body close to regular precession in the lagrange case. Izvestiya Akademii Nauk SSSR, Mekhanika Tverdogo Tela, 22(6):8–17, 1987. [21] D. D. Leshchenko and S. N. Sallam. Perturbed rotational motions of a rigid body sim- ilar to regular precession. Journal of Applied Mathematics and Mechanics, 54(2):183– 190, 1990. [22] T. S. Amer. On the rotational motion of a gyrostat with mass distribution. Nonlinear Dynamics, 54:189–198, 2008. [23] T. S. Amer and I. M. Abady. On the motion of a gyro in the presence of a newtonian T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 27 of 35 force field and applied moments. Mathematics and Mechanics of Solids, 23(9):1263– 1273, 2018. [24] A. I. Ismail, T. S. Amer, and S. A. El Banna. Electromagnetic gyroscopic motion. Journal of Applied Mathematics, 2012. [25] T. S. Amer. The rotational motion of the electromagnetic symmetric rigid body. Applied Mathematics and Information Sciences, 10(4):1453–1464, 2016. [26] V. S. Elfimov. Existence of periodic solutions of equations of motion of a solid body similar to the lagrange gyroscope. Journal of Applied Mathematics and Mechanics, 42(2):251–258, 1978. [27] A. I. Ismail. Applying the large parameter technique for solving a slow rotary motion of a disc about a fixed point. International Journal of Aerospace Engineering, 2020:7, 2020. [28] A. I. Ismail. Solving a problem of rotary motion for a heavy solid using the large parameter method. Advances in Astronomy, 2020:7, 2020. [29] E. T. Whittaker. A Treatise on the Analytical Dynamics of Particles and Rigid Bodies: With an Introduction to the Problem of Three Bodies. Cambridge University Press, 1988. [30] P. C. Hughes. Spacecraft Attitude Dynamics. Dover Publications, 2004. [31] J. J. Craig. Introduction to Robotics: Mechanics and Control. Springer, 2004. [32] T. S. Amer, W. S. Amer, M. Fakharany, A. H. Elneklawy, and H. F. El-Kafly. Mod- eling of the euler-poisson equations for rigid bodies in the context of the gyrostatic influences: An innovative methodology. European Journal of Pure and Applied Math- ematics, 18(1):5712, 2025. [33] T. S. Amer, H. A. El-Kafly, A. H. Elneklawy, and W. S. Amer. Modeling analysis on the influence of the gyrostatic moment on the motion of a charged rigid body subjected to constant axial torque. Journal of Low Frequency Noise, Vibration and Active Control, 43(4):1593–1610, 2024. [34] T. S. Amer, A. H. Elneklawy, and H. F. El-Kafly. Analysis of euler’s equations for a symmetric rigid body subject to time-dependent gyrostatic torque. Journal of Low Frequency Noise, Vibration and Active Control, 2025. [35] T. S. Amer, H. F. El-Kafly, A. H. Elneklawy, and A. A. Galal. Analyzing the dy- namics of a charged rotating rigid body under constant torques. Scientific Reports, 14(1):9839, 2024. [36] A. A. Galal, T. S. Amer, A. H. Elneklawy, and H. F. El-Kafly. Studying the influence of a gyrostatic moment on the motion of a charged rigid body containing a viscous incompressible liquid. The European Physical Journal Plus, 138:959, 2023. [37] F. M. El-Sabaa, T. S. Amer, H. M. Gad, and M. A. Bek. On the motion of a damped rigid body near resonances under the influence of harmonically external force and moments. Results in Physics, 19:103352, 2020. [38] J. Awrejcewicz, R. Starosta, and G. S. Kamińska. Asymptotic analysis of resonances in nonlinear vibrations of the 3-dof pendulum. Differential Equations and Dynamical Systems, 21:123–140, 2013. [39] T. S. Amer, M. A. Bek, and M. K. Abouhmr. On the vibrational analysis for the mo- T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 28 of 35 tion of a harmonically damped rigid body pendulum. Nonlinear Dynamics, 91:2485– 2502, 2018. [40] T. S. Amer. The dynamical behavior of a rigid body relative equilibrium position. Advances in Mathematical Physics, 2017:1–13, 2017. [41] R. W. Hamming. Numerical Methods for Scientists and Engineers. Dover Publica- tions, 1987. [42] I. M. Abady, T. S. Amer, H. M. Gad, and M. A. Bek. The asymptotic analysis and stability of 3dof non-linear damped rigid body pendulum near resonance. Ain Shams Engineering Journal, 91(2):101554, 2022. [43] T. S. Amer, F. M. El-Sabaa, S. K. Zakria, and A. A. Galal. The stability of 3-dof triple-rigid-body pendulum system near resonances. Nonlinear Dynamics, 110:1339– 1371, 2022. [44] S. H. Strogatz. Nonlinear Dynamics and Chaos: With Applications to Physics, Biol- ogy, Chemistry, and Engineering. Princeton University Press, Princeton, 2015. [45] R. Jarrar, R. Safdar, N. F. Mohammad, O. Florea, and J. Asad. Dynamical node of case study on mass spring system on a massless cart: Compared analytical and numerical solutions. Journal of Mechanics of Continua and Mathematical Sciences, 20:1–10, 2025. [46] O. A. Florea, D. Shahroor, R. Wannan, and J. Asad. Pendulum between two springs using ms-dtm. Nonlinear Dynamics, 100(2):025240, 2025. [47] T. Al-khader, D. K. Maiti, T. Roy, O. Florea, and J. Asad. Time dependent harmonic oscillator via om-hpm. Journal of Computational and Applied Mechanics, 56(1):264– 275, 2025. [48] T. Roy, A. Soqi, D. K. Maiti, R. Wannan, and J. Asad. Pendulum attached to a vibrating point: Semi-analytical solution by optimal and modified homotopy pertur- bation method. Alexandria Engineering Journal, 111:396–403, 2025. [49] T. S. Amer and I. M. Abady. On the application of kbm method for the 3-d motion of asymmetric rigid body. Nonlinear Dynamics, 89:1591–1609, 2017. [50] A. H. Nayfeh. Introduction to Perturbation Techniques. Wiley-VCH, 2011. [51] T. S. Amer, A. H. Elneklawy, and H. F. El-Kafly. Dynamical motion of a space- craft containing a slug and influenced by a gyrostatic moment and constant torques. Journal of Low Frequency Noise, Vibration and Active Control, 2025. [52] T. S. Amer and I. M. Abady. On the solutions of the euler’s dynamic equations for the motion of a rigid body. Journal of Aerospace Engineering, 30(4):04017021, 2017. [53] D. F. Griffiths and D. J. Higham. Numerical Methods for Ordinary Differential Equa- tions: Initial Value Problems. Springer, 2010. [54] T. S. Amer. Motion of a rigid body analogous to the case of euler and poinsot. Analysis, 24:305–315, 2004. [55] T. S. Amer, A. H. Elneklawy, and H. F. El-Kafly. A novel approach to solving euler’s nonlinear equations for a 3dof dynamical motion of a rigid body under gyrostatic and constant torques. Journal of Low Frequency Noise, Vibration and Active Control, 44(1):111–129, 2025. T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 29 of 35 Appendix 1 Γ11 = −1 2d 2 2[χ− kδ(1 + χδ)]− 1 br10 {A1[b(1 + χδ)− 1][d1 − δ + δd2(k + 1)] −A3[b(1 + χδ)− 1][d1χ− 1− χδ + d2(k + 1)(1 + χδ)] −A1A3[b(1 + χδ)− 1][χ+ δ(1 + χδ)(k + 1)]}+ λ́2 ρ1ρ2 {A2 b [b(1 + χδ)− 1] ×[ρ2(1 + χδ) + χδρ1 + δ2(1 + χδ)(ρ1 + ρ2)(k + 1)]− 1 2ρ1ρ2 δ2χλ́2r 2 10(1 +χδ)2(χ− (k + 1)(1 + χδ))(ρ1 + ρ2) 2}, Γ12 = −1 2δ 2[χ− kδ(1 + χδ)], Γ13 = −1 2(1 + χδ)2[χ− k(1 + χδ)], Γ14 = 1 br10 [b(1 + χδ))− 1][d1 − δ + δd2(k + 1)]− δd2[χ− kδ(1 + χδ)], Γ15 = 1 br10 [b(1 + χδ))− 1][d1χ− 1− χδ + d2(1 + k)(1 + χδ)]− d2(1 + χδ)[χ− kδ(1 + χδ)], Γ16 = 1 br10 [b(1 + χδ))− 1][χ+ δ(1 + k)(1 + χδ)]− δ(1 + χδ)[χ− kδ(1 + χδ)], Γ17 = λ́2 ρ1ρ2 { [b(1+χδ))−1] b [ρ2(1 + χδ) + χδρ1 + δ2(1 + χδ)(1 + k)(ρ1 + ρ2)] −δ2r10(1 + χδ)(ρ1 + ρ2)(χ− kδ(1 + χδ))}, Γ18 = λ́2 ρ1ρ2 { [b(1+χδ))−1] b {χ[ρ2(1 + χδ) + χδρ1] + δ(1 + k)(1 + χδ)2(ρ1 + ρ2)} −δλ́2r10(1 + χδ)2(ρ1 + ρ2)(χ− δk(1 + χδ))}, Γ21 = 1 2d 2 2[k(1 + χδ)2 − χ2]− 1 r10 {A1[χ(1 + χδ)][d1 − δ + δd2(1 + k)] −A3[χ(1 + χδ)− 1][d1χ− 1− χδ + d2(1 + k)(1 + χδ)]−A1A3[χ(1 + χδ)][χ+ δ(1 + χδ)(k + 1)]} + χλ́2 ρ1ρ2 {A2(1 + χδ)[ρ2(1 + χδ) + χδρ1 + δ2(1 + χδ)(1 + k)(ρ1 + ρ2)]− 1 2ρ1ρ2 δ2λ́22r 2 10(1 + χδ)2(ρ1 +ρ2) 2)(k(1 + χδ)2 − χ2)}, Γ22 = −1 2 [k(1 + χδ)2 − χ2]δ2, Γ23 = −1 2(1 + χδ)2[k(1 + χδ)2 − χ2], Γ24 = χ r10 {(1 + χδ)[d1 − δ + δd2(1 + k)] + δr10d2[k(1 + χδ)2 − χ2]}, Γ25 = χ(1+χδ) r10 {[d1χ+ (1 + χδ)(d2(1 + k)− 1)] + d2r10[k(1 + χδ)2 − χ2}, Γ26 = χ r10 (1 + χδ){[χ+ δ(1 + χδ)(1 + k)]− δr10[(1 + χδ)2 − χ2̂]}, Γ27 = χ(1+χδ)λ́2 ρ1ρ2 {[ρ2(1 + χδ) + χδρ1 + δ2(1 + χδ)(ρ1 + ρ2)(1 + k)] + δ2 χ [k(1 + χδ)2] −χ2]r10(ρ1 + ρ2)}, Γ28 = λ́2χ(1+χδ) ρ1ρ2 {[χ(ρ2(1 + χδ) + χδρ1) + δ(1 + k)(1 + χδ)2(ρ1 + ρ2)] + δr10(1 + χδ)2 ×[k(1 + χδ)2 − χ2](ρ1 + ρ2)}, T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 30 of 35 Γ31 = 1 br10 {A1(1 + χδ)λ́2[d1 − δ + δd2(1 + k)] +A3[(1 + χδ)λ́2[d1χ −1− χδ + d2(1 + χδ)(1 + k)]] +A1A3[(1 + χδ)λ́2[χ+ δ(1 + χδ)(1 + k)]]} − A2 bρ1ρ2 ×[(1 + χδ)λ́2 2 [ρ2(1 + χδ) + χδρ1 + δ2(1 + k)(1 + χδ)(ρ1 + ρ2)]], Γ32 = 0, Γ33 = 0, Γ34 = − λ́2 br10 (1 + χδ)[d1 + δ(d2(1 + k)− 1)], Γ35 = − λ́2 br10 {(1 + χδ)[d1χ+ (1 + χδ)(d2(1 + k)− 1)]}, Γ36 = − λ́2 br10 (1 + χδ)[χ+ δ(1 + k)(1 + χδ)], Γ37 = − λ́2 bρ1ρ2 (1 + χδ){(1 + χδ)[ρ2 + δ2(ρ1 + ρ2)(1 + k)] + χδρ1}, Γ38 = − λ́2 bρ1ρ2 (1 + χδ){χ[ρ2(1 + χδ) + χδρ1] + δ(1 + k)(1 + χδ)2(ρ1 + ρ2)}, Γ41 = 1 2d 2 2[χd1 + (1 + χδ)(1− kd2)]− 1 br10 {A1[ad1(1 + χδ)− χd2] ×[d1 − δ + δd2(k + 1)]−A3[ad1(1 + χδ)− χd2][d1χ− 1− χδ + d2(1 + k)(1 + χδ)] −A1A3[ad1(1 + χδ)− χd2][χ+ δ(1 + χδ)(1 + k)]}+ λ́2 ρ1ρ2 {A2 b [ad1(1 +χδ)− χd2][ρ2(1 + χδ) + χδρ1 + δ2(1 + χδ)(ρ1 + ρ2)(1 + k)] + 1 2ρ1ρ2 δ2λ́21r 2 10(1 + χδ)2 ×[χd1 + (1 + χδ)(1− kd2)](ρ1 + ρ2) 2}, Γ42 = 1 2 [χd1 + (1 + χδ)(1− kd2)]δ 2, Γ43 = 1 2 [χd1 + (1 + χδ)(1− kd2)](1 + χδ)2, Γ44 = 1 br10 {[ad1(1 + χδ)− χd2][d1 + δ(d2(k + 1)− 1)]}+ δ[χd1 + (1 + χδ)(1− kd2)]d2, Γ45 = 1 br10 {[ad1(1 + χδ)− χd2][d1χ+ (1 + χδ)(d2(k + 1)− 1)]}+ (1 + χδ)[χd1 +(1 + χδ)(1− kd2)]d2, Γ46 = 1 br10 {[ad1(1 + χδ)− χd2][χ+ δ(1 + χδ)(1 + k)]}+ δ(1 + χδ)[χd1 + (1 + χδ)(1− kd2)], Γ47 = λ́2 bρ1ρ2 {[ad1(1 + χδ)− χd2][(1 + χδ)(ρ2 + δ2(1 + k)(ρ1 + ρ2)) + χδρ1] +δ2(1 + k)br10(1 + χδ)(ρ1 + ρ2) + [d1χ+ (1 + χδ)(1− k)d2]}, Γ48 = λ́2 bρ1ρ2 {[ad1(1 + χδ)− χd2]{χ[ρ2(1 + χδ) + δχρ1] + δ(1 + χδ)2(ρ1 + ρ2(1 + k)) +δbr10(1 + χδ)2(ρ1 + ρ2)[χd1 + (1 + χδ)(1− kd2)]}, Γ51 = 1 2d 2 2(1− δ2k) + δ r10 {A1[(d1 + δ((1 + k)d2 − 1)] +A3[(d1χ− 1− χδ) + (1 + k)(1 + χδ)d2] +A1A3[χ+ (1 + χδ)δ(1 + k)]} − λ́2δ ρ1ρ2 {A2[ρ2(1 + χδ) + χδρ1 + δ2(1 + χδ)(ρ1 + ρ2)(k + 1)] + 1 2ρ1ρ2 δr210λ́ 2 2(1 + χδ)2(ρ1 + ρ2) 2(1− kδ2)}, Γ52 = 1 2δ 2(1− kδ2), T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 31 of 35 Γ53 = 1 2(1 + χδ)2(1− kδ2), Γ54 = δ r10 [d2(r10(1− kδ2)− δ(k + 1)) + δ − d1], Γ55 = d2(1 + χδ)(1− kδ2)− δ r10 [(1 + χδ)(d2(1 + k) + 1)− d1χ], Γ56 = δ r10 [(1 + χδ)(1− δ(1 + k))− χ] + δ(1 + χδ)(1− kδ2), Γ57 = λ́2δ ρ1ρ2 {[ρ2(1 + χδ) + χρ1 + δ(ρ1 + ρ2)(1 + χδ)(k + 1)]] + δr10(1 + χδ)(ρ1 + ρ2)(1− kδ2)}, Γ58 = λ́2δ ρ1ρ2 (ρ1 + ρ2){δr10(1− kδ2)(1 + χδ)2 − [(1 + χδ)2(1 + k) + χδ[ρ2(1 + χδ) + χδρ1]]}, Γ61 = 1 br10 {[d1 + δ(d2(k + 1)− 1)](1 + bχδ)A1 + (1 + bχδ)[(1 + χδ)(d2(k + 1)− 1) + d1χ]A3 +(1 + bχδ)[χ+ δ(k + 1)(1 + χδ)]A1A3} + d22 2 (χ− kδ(1 + δχ))− λ́2 bρ1ρ2 (1 + bχδ)[ρ2(1 + χδ) + χδρ1 + (ρ1 + ρ2) ×δ2(1 + χδ)(k + 1)]A3 + δ2 2 (χ− kδ(1 + χδ))( λ́2A2r10 ρ1ρ2 )2(1 + χδ)2(ρ1 + ρ2) 2], Γ62 = δ2 2 (χ− kδ(1 + χδ)), Γ63 = (1+χδ)2 2 (χ− kδ(1 + χδ)), Γ64 = − 1 br10 [d1 − δ(1− d2(1 + k))](1 + bχδ) + δ(χ− kδ(1 + χδ))d2, Γ65 = − 1 br10 (1 + bχδ)[(1 + χδ)(d2(k + 1)− 1) + χd1] + (1 + χδ)(χ− kδ(1 + χδ)))d2, Γ66 = − 1 br10 [χ+ δ(k + 1)(1 + χδ)](1 + bχδ) + δ(1 + χδ)(χ− kδ(1 + χδ)), Γ67 = − λ́2 bρ1ρ2 {(1 + bχδ)[ρ2(1 + χδ) + χδρ1 + δ2(1 + χδ)(ρ1 + ρ2)(k + 1)] −δ2br10(1 + χδ)(χ− kδ(1 + χδ))(ρ1 + ρ2)}, Γ68 = − λ́2 bρ1ρ2 {χ(1 + bχδ)[(ρ1 + ρ2)χδ + ρ2] + δ(1 + χδ)2(k + 1)(ρ1 + ρ2)] +δbr10(1 + χδ)2(χ− kδ(1 + χδ))(ρ1 + ρ2)}, Γ71 = − δλ́2[d1−δ+δd2(k+1)] br10 A1 − δλ́2[d1χ−(1+χδ)(d2(1+k)−1)] br10 A3 − δλ́2[χ+δ(1+k)(1+χδ)] br10 A1A3 + δλ́2 2 bρ1ρ2 [δχρ1 + (1 + χδ)ρ2 + δ2 ×(1 + χδ)(ρ1 + ρ2)(k + 1)]A2, Γ72 = 0, Γ73 = 0, Γ74 = δλ́2[(d1−δ)+δd2(k+1)] br10 , Γ75 = δλ́2[d2(1+χδ)(k+1)+(d1χ−(1+χδ))] br10 , Γ76 = δλ́2 br10 [χ+ δ(1 + χδ)(k + 1)], Γ77 = δλ́2 2 bρ1ρ2 [χδρ1 + ρ2(1 + χδ) + δ2(1 + χδ)(k + 1)(ρ1 + ρ2)], Γ78 = δλ́2 2 bρ1ρ2 {χ[χδρ1 + ρ2(1 + χδ)] + δ(1 + χδ)2(ρ1 + ρ2)(k + 1)}, T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 32 of 35 Γ81 = − [(d1−δ)+δd2(k+1)][d2−δad1)] br10 A1 − 1 br10 [(d1χ− (1 + χδ)) +d2(1 + χδ)(k + 1)][d2 − δad1]A3 − 1 2(d1 + δ(1− kd2))d 2 2 + [χ+δ(1+χδ)(k+1)][d2−δad1] br10 A1A3 + λ́2 bρ1ρ2 [χδρ1 + ρ2(1 + χδ) + δ2 ×(1 + χδ)(ρ1 + ρ2)(k + 1)][d2 − δad1]A2 − 1 2 δ2λ́2 2 r210(1+χδ)2 (ρ1ρ2)2 (d1 +δ(1− kd2))(ρ1 + ρ2) 2, Γ82 = − δ2 2 (d1 + δ(1− kd2)), Γ83 = − (1+χδ)2 2 (d1 + δ(1− kd2)), Γ84 = [(d1−δ)+δd2][d2−δ(ad1+λ́1)] br10 − δd2(d1 + δ), Γ85 = 1 br10 [(d1χ− (1 + χδ)) + (1 + χδ)d2(k + 1)][d2 − δad1]− d2(1 + χδ) ×(d1 + δ(1− kd2)), Γ86 = [χ+δ(k+1)(1+χδ)][d2−δad1] br10 − δ(1 + χδ)(d1 + δ(1− kd2)), Γ87 = λ́2 bρ1ρ2 [(1 + χδ)ρ2 + χδρ1 + δ2(k + 1)(1 + χδ)(ρ1 + ρ2)][d2 − δad1] − δ2λ́2r10 ρ1ρ2 (1 + χδ)(ρ1 + ρ2)(d1 + δ(1− kd2)), Γ88 = λ́2 bρ1ρ2 {χ[ρ1χδ + (1 + χδ)ρ2] + δ(1 + χδ)2(ρ1 + ρ2)(k + 1)}[d2 − δad1] − δλ́2r10(1+χδ)2 ρ1ρ2 (ρ1 + ρ2)(d1 + δ(1− kd2)), T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 33 of 35 Appendix2 G10 = Γ41 + 1 2 [(2Γ42 − Γ41)(A 2 1 +A2 2) + (2Γ43 − Γ25)A 2 3]− 2Γ17A1A2, G11 = Γ44 − 2Γ11 − 2(A2 1 +A2 2)Γ12 − 1 2 [4Γ13 + 1 ρ1 (ρ1 + ρ2)Γ26]A 2 3 − 2 ×ρ1µ0 − λ́2r10(1+χδ) ρ1 Γ17 − λ́2r10δ ρ2 Γ27 + Γ37, G12 = 1 2ρ1 (ρ1 − ρ2)A 2 3Γ26, G13 = −3 2(A 2 1 −A2 2)Γ14 +A1A2Γ17, G14 = −3A1A2Γ14 − 1 2Γ17(3A 2 1 −A2 2)− 1 2ρ1 (ρ1 + ρ2)A1A2Γ27, G15 = Γ45 − 1 2(3− ρ2 ρ1 )(A2 1 −A2 2)Γ16 − 1 ρ1 (ρ1 + ρ2)[Γ21 +A2 3Γ23 +(A2 1 +A2 2)Γ22]− λ́2r10δ ρ1 Γ28 + ρ2 ρ1 [Γ38 − λ́2r10(1+χδ) ρ1 Γ18], G16 = Γ46 − 1 2 [2Γ15 + 1 ρ1 (ρ1 + ρ2)Γ24], G17 = −1 2(Γ24 + ρ2 ρ1 Γ15), G18 = −1 2 [2Γ15 + 1 ρ1 (ρ1 + ρ2)Γ24), G19 = ρ2 2ρ1 Γ15 + 1 2Γ24, G110 = 1 2Γ61A3[ 1 ρ1 (ρ1 − ρ2)A 2 1 − 2(A2 1 +A2 2)], G111 = −(3− ρ2 ρ1 )A1A2A3Γ16, G112 = − 1 2ρ1 (ρ1 + ρ2)Γ26, G113 = − 1 2ρ1 (ρ1 − ρ2)Γ26, G114 = − 1 2ρ1 (ρ1 + ρ2)A 2 3Γ25, G115 = − 1 2ρ1 (ρ1 + 2ρ2)A 2 3Γ28, G116 = − 1 ρ1 (ρ1 − ρ2)Γ36 − (1 + ρ2 2ρ1 )Γ18 + 1 2ρ1 (ρ1 − ρ2) λ́2r10δ ρ2 Γ26, G117 = −1 2 [ 1 ρ1 (ρ1 − ρ2)− 1]Γ27, G118 = −1 2Γ27 − 1 2ρ1 (ρ1 + ρ2)Γ27 − Γ18, G119 = − ρ2 2ρ1 Γ18, G120 = − 1 2ρ1 (ρ1 + ρ2)Γ27, G121 = Γ47 + λ́2r10δ ρ2 Γ24 + λ́2r10(1+χδ) ρ1 Γ14 − Γ34, G122 = Γ48 + λ́2r10δ ρ1 Γ25 + ρ2 ρ1 λ́2r10(1+χδ) ρ1 Γ15 − ρ2 ρ1 Γ35, G123 = (1− ρ2 ρ1 ) λ́2r10(1+χδ) ρ1 A1A3Γ16, T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 34 of 35 G20 = −2r10λ́2µ0(1 + χδ) + [Γ11 + (A2 1 +A2 2)Γ12 +A2 3Γ13] λ́2r10(1+χδ) ρ1 + 1 2A 2 3Γ28 + [Γ21 +(A2 1 +A2 2)Γ22 +A2 3Γ23] λ́2r10δ ρ2 − [Γ31 + (A2 1 +A2 2)Γ32 +A2 3Γ33]− (A2 1 +A2 2)Γ17, G21 = −2ρ1µ0 − 2[Γ11 + (A2 1 +A2 2)Γ12 +A2 3Γ13]− λ́2r10δ ρ2 Γ27 − 1 2(1− ρ2 ρ1 )A2 3Γ26 +Γ37 + Γ44 − 2λ́2r10(1+χδ) ρ1 Γ17, G22 = −Γ34 + Γ47 + λ́2r10δ ρ2 Γ24 − 2(1+χδ)λ́2r10 ρ1 Γ14, G23 = 1 2ρ1 (ρ1 + ρ2)A 2 3Γ26, G24 = ρ2 2ρ1 (A2 1 +A2 2)Γ16 + ρ2 ρ1 Γ45 − λ́2r10δ ρ2 Γ28 + 2(1+χδ)λ́2r10 ρ1 Γ18 + Γ38 −(1 + ρ2 ρ1 )[Γ21 + (A2 1 +A2 2)Γ22 +A2 3Γ23], G25 = 1 2(Γ15 − Γ24) + 1 ρ1 (ρ1 − ρ2)Γ46, G26 = 1 2 [Γ18 − (1 + ρ2 ρ1 )Γ27], G27 = −1 2(Γ27 + ρ2 ρ1 )Γ18) + λ́2r10δ ρ2 Γ26 − 2(1+χδ)λ́2r10 ρ1 Γ16 − Γ36, G28 = −1 2(1 + ρ2 ρ1 )Γ24 − ρ2 2ρ1 Γ15, G29 = 1 2(Γ15 − Γ24), G210 = 1 2ρ1 [(ρ1 + ρ2)Γ24 + ρ2Γ15], G211 = 1 2 [−Γ18 + (1 + ρ2 ρ1 )Γ27], G212 = −1 2 [ ρ2 ρ1 Γ18 + Γ27], G213 = 1 2ρ1 (2ρ1 − ρ2)(A 2 2 −A2 1)A3Γ16, G214 = − 1 2ρ1 (ρ1 + ρ2)Γ26, G215 = 1 2ρ1 (ρ1 − ρ2)Γ26, G216 = 1 2A 2 3Γ25, G217 = 1 ρ1 (ρ1 − ρ2)A1A2Γ16 − ρ2 ρ1 Γ18 + δλ́2r10 ρ2 Γ25 − 2(1+χδ)λ́2r10 ρ1 Γ15 − Γ35, G218 = −1 2(1 + 2ρ2 ρ1 )Γ28, G32 = Γ85 − 2(ρ2µ0 + Γ61 + Γ63A 2 3)− δλ́2r10 ρ2 Γ68 − (1+χδ)λ́2r10 ρ1 Γ58 + Γ78, G33 = −3 2Γ68, G34 = −1 2Γ65, G40 = 2δλ́2r10µ0 + δλ́2r10 ρ2 [Γ61 + (A2 1 +A2 2)Γ62 +A2 3Γ63] + (1 + ρ1 ρ2 ) (1+χδ)λ́2r10 ρ1 − (1+χδ)λ́2r10 ρ2 [Γ51 + (A2 1 +A2 2)Γ52 +A2 3Γ53]− [Γ71 + (A2 1 +A2 2)Γ72 +A2 3Γ73] +1 2A 2 3Γ68]− ρ1 2ρ2 (A2 1 +A2 2)Γ57, G41 = − δλ́2r10 ρ2 Γ67 − (1 + ρ1 ρ2 ) + (1+χδ)λ́2r10 ρ2 Γ57 + 1 2(1− ρ1 ρ2 )A2 3Γ66 +Γ77 + ρ1 ρ2 Γ84, G42 = A2 3Γ66, G43 = δλ́2r10 ρ2 Γ64 − (1+χδ)λ́2r10 ρ2 Γ54 + ρ1 ρ2 Γ87 − Γ74, G44 = − ρ1 2ρ2 (A2 1 −A2 2)Γ54 + ρ1 ρ2 A1A2Γ57, G45 = − ρ1 2ρ2 (A2 1 −A2 2)Γ57 − ρ1 ρ2 A1A2Γ54, G46 = −4ρ2µ0 − 2[Γ61 + (A2 1 +A2 2)Γ62 +A2 3Γ63] + δλ́2r10 ρ2 (Γ56 −2Γ68) + Γ85 + Γ78 − 1 2(1− ρ1 ρ2 )(A2 1 +A2 2)Γ56 + (1+χδ)λ́2r10 ρ1 Γ58, G47 = − ρ1 2ρ2 Γ64 − (1− ρ1 ρ2 )Γ86, G48 = −1 2(2Γ64 + Γ55), G49 = δλ́2r10 ρ2 Γ66 − (Γ76 + 1 2Γ58)− (1+χδ)λ́2r10 ρ2 Γ56 − ρ1 2ρ2 Γ67, T. S. Amer et al. / Eur. J. Pure Appl. Math, 18 (3) (2025), 6281 35 of 35 G410 = −Γ67, G411 = − ρ1 2ρ2 Γ64, G412 = 1 2(2Γ64 + Γ55), G413 = −1 2(Γ58 + ρ1 ρ2 Γ67), G414 = −Γ67, G415 = − 1 2ρ2 (ρ2 − ρ1)(A 2 1 −A2 2)A3Γ56, G416 = 0, G417 = 1 2ρ2 (ρ2 − ρ1)Γ66, G418 = −Γ66, G419 = −3 2A 2 3Γ65, G420 = δλ́2r10 ρ2 Γ65 + Γ88 − Γ75 − (1+χδ)λ́2r10 ρ1 Γ55 + (1− ρ1 ρ2 )A1A2Γ56, G421 = −3 2A 2 3Γ68,