[BACK] Sound & Vibration DOI: 10.32604/sv.2022.014166

ARTICLE

Damped Mathieu Equation with a Modulation Property of the Homotopy Perturbation Method

Department of Mathematics, Faculty of Education, Ain Shams University, Roxy, Cairo, Egypt
*Corresponding Author: Nasser S. Elgazery. Email: nasser522000@gmail.com
Received: 06 September 2020; Accepted: 25 June 2021

Abstract: In this article, the main objective is to employ the homotopy perturbation method (HPM) as an alternative to classical perturbation methods for solving nonlinear equations having periodic coefficients. As a simple example, the nonlinear damping Mathieu equation has been investigated. In this investigation, two nonlinear solvability conditions are imposed. One of them was imposed in the first-order homotopy perturbation and used to study the stability behavior at resonance and non-resonance cases. The next level of the perturbation approaches another solvability condition and is applied to obtain the unknowns become clear in the solution for the first-order solvability condition. The approach assumed here is so significant for solving many parametric nonlinear equations that arise within the engineering and nonlinear science.

Keywords: Damped Mathieu Equation; parametric nonlinear oscillator; resonance instability; homotopy perturbation method (HPM)

1  Introduction

In wide engineering and physical applications, the nonlinear oscillators exist. Also, the parametric excitation takes place when a modifying physical parameter, such as a moment of stiffness or inertia, acts in a forcing model. This excitement yields a variable time coefficient, commonly an oscillation, in the governing system of motion. On the other hand, an external excitation outcome acting as an inhomogeneous part in the model of motion. Furthermore, minor parametric excitement produces a major response when the frequency of the excitement is far from the fundamental resonance, as shown in .

A classical example of parametric excitation is the swinging pendulum with oscillating support. The equation of motion describing the model is the well-known Mathieu equation. In 1868 Mathieu studied the vibration of elliptical membranes . Consequently, he introduced the Mathieu equation that is an example of a linear differential equation (LDE) with parametric excitation. The Mathieu equation has application to the dynamics of passive towed arrays in submarines, as well as serving as a useful model for many interesting problems in physics, biology, applied mathematics, and engineering mechanics fields . For some of the non-linear variations of the Mathieu, the equation has been presented in [8,9]. Moreover, the oscillations of the mechanical systems under the action of an oscillatory external force may reveal a Duffing problem, for instance, see references . Recently, Moatimid  attempted to study the stability analysis of a parametric Duffing oscillator. In this investigation Moatimid showed that the cubic stiffness parameter and the damped parameter have a destabilizing influence, however, the parametric and natural frequencies are of stabilizing influences.

The main target in the present work is how to achieve accurate approximate solutions of the nonlinear oscillators with highly strong nonlinearity. In recent centuries, many analytical approaches were developed to work out the periodic motion of nonlinear oscillators, such as the averaging method, perturbation methods, harmonic balance method, and the generalized harmonic method. The classical perturbation procedure depends on small parameters and chooses unsuitable small parameters that can lead to wrong solutions . Therefore, a new perturbation technique was first proposed by He et al. . This technique is named as the HPM, which represents a combination of the Homotopy analysis and classical perturbation methods. It has a full promise of the traditional perturbation techniques. The major property of the HPM is in its ability and flexibility to deal with many types of linear and nonlinear differential equations conveniently and accurately. Further, the HPM provides us with an appropriate direction to calculate an approximate or an analytic solution to several models arising in different fields. He  was built the most two considerable steps in the criteria of the HPM with a suitable initial guess and suggested an alternative approach to the construction of the Homotopy equation. Hence, He applied HPM to solve the Lighthill equation , Duffing equation , and Blasius equation , then the idea goes through and has been applied to solve nonlinear wave equations , boundary value problems . Babolian et al.  applied the homotopy perturbation method to solve the Burgers, the modified Korteweg-de Vries, and regularized long-wave equations.

On the other hand, the HPM has more improved and developed by many engineers and scientists, for instance, a couple of the Laplace transforms and Homotopy perturbation method was implemented by El-Dib et al. . The HPM with two expanding parameters that efficient for some partial nonlinear equations was suggested by He  and El-Dib . Also, El-Dib  introduced a modified version of the HPM via the multiple scales technique. This new modification works particularly well for the nonlinear oscillators. Furthermore, away from the traditional HPM, Ren et al.  made another couple of the HPM and multiple time scales to become a powerful mathematical tool for many nonlinear equations. They displayed that the present procedure may be further afflicted by incorporating several known technologies. It provides solutions to nonlinear equations, whilst the classical perturbation technique became unsuccessful. Moreover, Rabbani  introduced a new homotopy perturbation approach for solving main non-linear models through the projection method. A new homotopy perturbation technique for solving linear and nonlinear Schrödinger equations has been addressed by Ayati et al. . Further, by utilizing the HPM, a novel approach in examining the nonlinear Rayleigh-Taylor instability is conducted by El-Dib et al. . Recently, a periodic solution of the cubic nonlinear Klein–Gordon equation using the He-multiple-scales method has been investigated by El-Dib . Also, El-Dib et al.  investigated the impact of fractional derivative properties on the periodic analytic solution of the nonlinear oscillations using the HPM. Moreover, He’s-multiple-scale scale to analyze the cubic-quintic Duffing equation has been analyzed by El-Dib et al. . El-Dib  presented a stability approach of a fractional-delayed Duffing oscillator. A Nonlinear Instability of a Cylindrical Interface between two MHD Darcian flows has been studied by Moatimid et al. . Further, El-Dib  introduced a modified multiple scale technique for the stability of the fractional delayed nonlinear oscillator. Besides, a periodic solution of the fractional sine-Gordon equation has been studied by Shen et al. . Elgazery  applied the HPM to give a periodic solution of the Newell-Whitehead-Segel model. Further, for more very useful modification of the homotopy perturbation approach, Yu et al.  introduced HPM with an auxiliary parameter for nonlinear oscillators. Also, HPM for Fangzhu oscillator has been used by He et al. . Finally, for more very useful modification of the homotopy perturbation approach, see .

Motivated by possibility applications in engineering, biology, and physics, which is based on studying the solution of the damped Mathieu equation. Hence, in the present work, our objective is to apply the HPM to linear or nonlinear equations having periodic coefficients such Mathieu equation which has been of great importance among researchers. For the presentation of this article; the rest of the manuscript is systematized as follows: Section 2 is introducing the HPM for the mathematical formulation. The modulation procedure, in detail, is displayed in Section 3. The non-resonance case, stability analysis of the non-resonance case, the resonance case of Ω is near ω , stability analysis of the linear Mathieu equation, and stability analysis for the nonlinear case are represented through Sections 4 to 8, respectively. Finally, the main obtained outcomes are summarized as concluding remarks in Section 9.

2  Mathematical Formulation

To explain the proposed technique, consider the following parametric pendulum equation as an illustrative example:

d2ydt2+μdydt+(ω2+2qcos2Ωt)siny=0, (1)

where μ,ω,q and Ω are real physical constants. μ is the damping coefficient parameter, ω2 is the square of the natural frequency of the unforced system, q is the amplitude of the parametric excitation, and Ω is the frequency of the parametric excitation. By using the Taylor expansion up to cubic power, the above parametric pendulum equation becomes the following damping cubic nonlinear Mathieu equation:

d2ydt2+μdydt+(ω2+2qcos2Ωt)y=(ω23!+2q3!cos2Ωt)y3+.... (2)

The homotopy perturbation method can be considered as a combination of the classical perturbation technique and the homotopy (whose origin is in the topology ), but not restricted to the limitations of traditional perturbation methods. For example, this method does require neither small parameter nor linearization and only requires little iteration to obtain accurate solutions  and .

We define the two parts of Eq. (2) as L(y) and N(y) , where

L(y)=d2ydt2+ω2y, and N(y)=μdydt+2qycos2Ωtky3, (3)

where, k=13!(ω2+2qcos2Ωt) .

Construct the homotopy statement as

H(y,ρ)=L(y)+ρN(y)=0;ρ[0,1]. (4)

As in He’s a homotopy perturbation method, it is obvious that when ρ=0 , Eq. (4) becomes the harmonic equation; L(y)=0. Thus,

d2y(t)dt2+ω2y(t)=0. (5)

According to linear differential equations theory, the general solution of Eq. (5) is expressed in terms of two linearly independent solutions, say, eiωt and eiωt . Thus, the composite solutions may be in the form

y(t)=Aeiωt+A¯eiωt, (6)

where A and its complex conjugate A¯ are arbitrary constants of integration. Eq. (4) becomes the original nonlinear Mathieu Eq. (2) as ρ=1 . For arbitrary the small parameter ρ , the solution of Eq. (4) can be sought in terms of ρ so that the function y(t) becomes y(t,ρ) . Accordingly, Eq. (4) can rewrite as

H(y,ρ)=(d2dt2+ω2)y(t,ρ)+ρ(μddt+2qcos2Ωtky2(t,ρ))y(t,ρ)=0. (7)

It can be noticed that the homotopy function (7) is essentially the same as (4), except for the function y(t,ρ) , which contains embedded the homotopy parameter ρ . The introduction of that parameter within the differential equation is a strategy to redistribute the periodic part between the successive iterations of the homotopy method, and thus increase the probabilities of finding the sought solution. Thus, as ρ moves from 0 to 1, the function y(t,ρ) moves from y0(t) to yapp(t) . Expand the function y(t,ρ) as a power series in the small parameter ρ such that

y(t,ρ)=y0(t)+ρy1(t)+ρ2y2(t)+..., (8)

where yn(t);n=1,2,3,... are unknowns in needs to determine. Substituting this expansion into the homotopy Eq. (7) and equating terms with identical powers of ρ , leads to

ρ0:d2y0dt2+ω2y0=0, (9)

ρ1:d2y1dt2+ω2y1=μdy0dt2qcos2Ωty0+ky03, (10)

ρ2:d2y2dt2+ω2y2=μdy1dt2qcos2Ωty1+3ky02y1,.... (11)

Eq. (9) can be satisfied by

y0(t)=Aeiωt+A¯eiωt. (12)

Substituting (12) into Eq. (10) gets the form

d2y1dt2+ω2y1=iωμ(AeiωtA¯eiωt)qA(ei(ω+2Ω)t+ei(ω2Ω)t)qA¯(ei(ω+2Ω)t+ei(ω2Ω)t)+k(A3e3iωt+3A2A¯eiωt+3A¯2Aeiωt+A¯3e3iωt). (13)

Before analyzing the first-order problem, we must distinguish between the two cases. The case of the frequency Ω doesn’t equal the nature frequency ω (which is known as the non-resonance case). The second one is the specific case when Ω approaches ω (which is known as the resonance case).

For arbitrary frequency Ω , there are secular terms that appear in the Eq. (13). Elimination such secular terms require that the arbitrary constant A be zero. This means that the expansion (8) cannot be successful to obtain a valid solution for exciting homotopy Eq. (7).

3  The Modulation Procedure

To obtain uniform expansions for problems of this kind, the expansion (8) needs to be modified. If we modulate the initial solution (6) so that the constant A becomes A(τ) with τ=ρt, such that

Then Eq. (12) in the modulate case becomes

Y0(t,τ)=A(τ)u0(t)+A¯(τ)u¯0(t), (15)

where

u0(t)=eiωt and u¯0(t)=eiωt. (16)

Consequently, the homotopy state, Eq. (7), in the modulated form becomes

(d2dt2+ω2)Y(t,τ,ρ)+ρ(μddt+2qcos2Ωt)Y(t,τ,ρ)=ρkY3(t,τ,ρ). (17)

It is convenient to choose the modulated function Y(t,τ,ρ) in separated variables as

Y(t,τ,ρ)=A(τ)u(t,ρ)+A¯(τ)u¯(t,ρ). (18)

The function u(t,ρ) can be expanded as a power series in the small parameter ρ such that

u(t,ρ)=u0(t)+ρu1(t)+ρ2u2(t)+..., (19)

where un(t);n=1,2,3,... are unknowns to be evaluated. If the expansion (19) is substituted into Eq. (18) then gets

Y(t,τ,ρ)=A(τ)(u0(t)+ρu1(t)+ρ2u2(t)+...)+cc.=Y0(t,τ)+ρY1(t,τ)+ρ2Y2(t,τ)+..., (20)

where cc. indicates to the complex conjugate for the preceding terms and

Yn(t,τ)=A(τ)un(t)+A¯(τ)u¯n(t). (21)

It is noted that:

ddtY(t,τ,ρ)=ddt[A(τ)u(t,ρ)+cc]=A(τ)u˙(t,ρ)+ρu(t,ρ)A(τ)+cc, (22)

and

d2dt2Y(t,τ,ρ)=A(τ)u¨(t,ρ)+2ρA(τ)u˙(t,ρ)+ρ2u(t,ρ)A(τ)+cc, (23)

where dots indicate differentiation concerning the time t , while dashes refer to the derivative for the time modulate τ. Substituting (18) into Eq. (17) using (22) and (23) gives

A(u¨+ω2u)+ρ(2u˙A+μAu˙+2qAucos2Ωt)+ρ2(A+μA)uρk(A3u3+3A2A¯u2u¯)+cc.=0. (24)

Eq. (24) remains to obey the same homotopy concept because it’s become the same harmonic Eq. (5) as ρ0 . limρ1A=ddτ(limρ1A)=0, besides, consequently the original Eq. (2) is found.

In the light of Eq. (19), the modulate homotopy Eq. (24) will be expanded as a power series in ρ so that the following non-homogenous harmonic equations are imposed

ρ0:A(u¨0+ω2u0)+A¯(u¨¯0+ω2u¯0)=0, (25)

ρ1:A(u¨1+ω2u1)+2u˙0A+μAu˙0+2qAu0cos2ΩtkA2(Au03+3A¯u02u¯0)+cc.=0, (26)

ρ2:A(u¨2+ω2u2)+2u˙1A+μAu˙1+2qAu1cos2Ωt+u0A+μAu03kA2[Au02u1+A¯(u02u¯1+2u0u¯0u1)]+cc.=0. (27)

It is noted that Eq. (25) has been satisfied by Eq. (16) and the zero-order solution for Eq. (17) as approved in Eq. (15). Substituting Eq. (16) into Eq. (26) becomes

A(u¨1+ω2u1)+[iω(2A+μA)3kA2A¯]eiωt+qA(ei(ω+2Ω)t+ei(ω2Ω)t)kA3e3iωt+cc.=0. (28)

This equation contains secular terms at the non-resonance case and other secular terms when the applied frequency Ω approaches the natural frequency ω .

4  The Non-Resonance Case

The analysis in this case concerned with the arbitrary chosen for the applied frequency Ω , in Eq. (28). At this stage, secular terms are removed when

A+12μA+3ik2ωA2A¯=0, (29)

with its complex conjugate one. This leads to obtaining the valid function u1(t) as

u1(t)=q4Ω(ei(ω+2Ω)t(Ω+ω)+ei(ω2Ω)t(Ωω))+kA28ω2e3iωt. (30)

Consequently, the solution of the first-order problem is formulated as

Y1(t,τ)=qA(τ)4Ω(ei(ω+2Ω)t(ω+Ω)+ei(ω2Ω)t(ωΩ))+qA¯(τ)4Ω(ei(ω+2Ω)t(ω+Ω)+ei(ω2Ω)t(ωΩ))+k8ω2(A3(τ)e3iωt+A¯3(τ)e3iωt). (31)

Substituting Eqs. (16) and (30) into Eq. (27), using Eq. (29), yields

A(u¨2+ω2u2)+[A+μA+q2A2(Ω2ω2)3k28ω2A3A¯2]eiωt+3k2A4A¯8ω2e3iωt3k2A58ω2e5iωt+q2A4Ω(Ω2ω2)[(Ωω)ei(ω+4Ω)t+(Ω+ω)ei(ω4Ω)t]+3kqA2A¯2ω(Ω2ω2)[(Ω2ω)ei(ω+2Ω)t(Ω+2ω)ei(ω2Ω)t]3kqA34Ω(Ω2ω2)[(Ωω)ei(3ω+2Ω)t+(Ω+ω)ei(3ω2Ω)t]+cc.=0. (32)

The valid solution requires to be removed the terms that producing unbounded solution. These terms imply the following nonlinear solvability condition:

A+μA+q2A2(Ω2ω2)3k28ω2A3A¯2=0. (33)

The second-order solution is found to be

Y2(t,τ)=q2A(τ)32Ω2(Ω2ω2)[(Ωω)(ω+2Ω)ei(ω+4Ω)t+(Ω+ω)(ω2Ω)ei(ω4Ω)t]3kqA3(τ)16Ω(Ω2ω2)2[(Ωω)2(2ω+Ω)ei(3ω+2Ω)t+(Ω+ω)2(Ω2ω)ei(3ω2Ω)t]+3kqA2(τ)A¯(τ)8Ω2ω(Ω2ω2)[(Ω2ω)(ω+Ω)ei(ω+2Ω)t+(Ω+2ω)(Ωω)ei(ω2Ω)t]k2A5(τ)64ω2e5iωt+3k2A4(τ)A¯(τ)64ω2e3iωt+cc. (34)

If the accuracy to the second-order perturbation is enough, then the approximate solution at the non-resonance case is formulated by substituting Eqs. (15), (16), (31) and (34) into Eq. (20), and setting ρ=1 , gets

Y(t)=limρ1,τtY(t,τ,ρ)=A(t)eiωt+k8ω2A3(t)e3iωtk2A5(t)64ω2e5iωt+3k2A4(t)A¯(t)64ω2e3iωt+qA(t)4Ω(ei(ω+2Ω)t(ω+Ω)+ei(ω2Ω)t(ωΩ))+3kqA2(t)A¯(t)8Ω2ω(Ω2ω2)[(Ω2ω)(ω+Ω)ei(ω+2Ω)t+(Ω+2ω)(Ωω)ei(ω2Ω)t]+q2A(t)32Ω2(Ω2ω2)[(Ωω)(ω+2Ω)ei(ω+4Ω)t+(Ω+ω)(ω2Ω)ei(ω4Ω)t]3kqA3(t)16Ω(Ω2ω2)2[(Ωω)2(2ω+Ω)ei(3ω+2Ω)t+(Ω+ω)2(Ω2ω)ei(3ω2Ω)t]+cc. (35)

5  Stability Analysis for the Non-Resonance Case

The stability criteria in the non-resonance case can be obtained from solving Eq. (29). One may use the following polar form :

A(τ)=12ξ(τ)eiη(τ), (36)

with real the unknown functions ξ(τ) and η(τ) . Insert Eq. (36) into the first-order solvability condition (29) which will separate into real and imaginary parts and gives

ξ(τ)=ξ0e12μτ and η(τ)=3k2μωξ0e12μτ+η0, (37)

where, ξ0 and η0 are integration constants. The stability criteria in the non-resonance case require that μ>0.

6  The Resonance Case Ω is Near ω

Return to the first-order problem Eq. (28) and re-analyzed it because of the nearness of Ω to ω . We express this approach by introducing the detuning parameter σ  such that

Ω=ω+ρσ. (38)

Accordingly, we have

i(ω2Ω)t=iωt+2iστ. (39)

Elimination of secular terms from Eq. (28), because of Eqs. (38) and (39) yields

A+12μAiq2ωA¯e2iστ+3ik2ωA2A¯=0. (40)

The first-order solution in this case is

Y1(t,τ)=q4Ω(ω+Ω)(Aei(ω+2Ω)t+A¯ei(ω+2Ω)t)+k8ω2(A3e3iωt+A¯3e3iωt). (41)

Using Eq. (41) with Eq. (27), we obtain the uniform solution for the second-order problem, and the following solvability is presented:

A+μA+q24Ω(ω+Ω)A+kq8ω2A3e2iστ3kq4Ω(ω+Ω)AA¯2e2iσt3k28ω2A3A¯2=0, (42)

with its complex conjugate. The valid function Y2(t,τ) is given by

Y2(t,τ)=3ik64ω3(2A+μA+2ikA2A¯)A2e3iωt3k2192ω4A5e5iωt+q232Ω2(ω+Ω)(ω+2Ω)Aei(ω+4Ω)t+iq8Ω2(ω+Ω)2[(A+12μA)(ω+2Ω)+3ikA2A¯]ei(ω+2Ω)tkq16(2ω+Ω)(ω+Ω)[3Ω(ω+Ω)12ω2]A3ei(3ω+2Ω)t+cc. (43)

The approximate solution up to the second-order is formulated by substituting from Eqs. (15), (16), (41) and (43) into Eq. (16) gets

Y(t)=limρ1τt(Y0+ρY1+ρ2Y2)=Aeiωt+A¯eiωt+k8ω2[A+3i8ω(2A+μA+2ikA2A¯)]A2e3iωt3k2192ω4A5e5iωt+q4Ω(ω+Ω){A+iq2Ω(ω+Ω)[(A+12μA)(ω+2Ω)+3ikA2A¯]}ei(ω+2Ω)t+q232Ω2(ω+Ω)(ω+2Ω)Aei(ω+4Ω)tkq16(2ω+Ω)(ω+Ω)[3Ω(ω+Ω)12ω2]A3ei(3ω+2Ω)t+cc. (44)

7  Stability Analysis of the Linear Mathieu Equation

In the limiting case as k0 into Eq. (2), linear damping Mathieu equation arrived. In this case, the two solvability conditions (40) and (42) that produced at the resonance case of Ω is near ω having the following limit case:

A+12μAiq2ωA¯e2iστ=0, (45)

A+μA+q24Ω(ω+Ω)A=0. (46)

The first-order solvability condition (45) can be used to find the stability picture at the resonance case. The second-order solvability condition (46) can be used to find the value of the detuning parameter σ .

It is easy to show that the Eq. (45) can be satisfied by the form

A(τ)=[(σ+q2ω)sinΘτ+iΘcosΘτ]e(iσ12μ)τ, (47)

where, the parameter μ must be positive, to find a damping solution. The argument Θ is given by the following characteristic equation:

Θ2=σ2q24ω2. (48)

The parameter σ can be evaluated by substituting Eq. (47) into the second-order solvability condition (46) to gets

σ=q2ωorσ=±ω4q(μ2q2ω2q2Ω(Ω+ω)). (49)

The use of the first value of σ , Eq. (48) yields a zero solution for Eq. (45). For a non-zero solution, the other values σ are conforming. Inserting Eq. (49) into Eq. (48) gets

Θ2=ω216q2(μ2q2ω2q2Ω(Ω+ω))2q24ω2. (50)

The stability criteria require that the right-hand-side of Eq. (50) be positive, which implies that

(μ2q2ω2q2Ω(Ω+ω))24q4ω4>0. (51)

Stability condition (51) can be rearranged in powers of the applied frequency Ω as

Ω2(μ2ω23q2)+Ωω(μ2ω23q2)q2ω2>0, (52)

and

Ω2(μ2ω2+q2)+Ωω(μ2ω2+q2)q2ω2<0. (53)

The transition curves separating stable state from an unstable state corresponding to

Ω1=ω(μ2ω23q2)ωμ4ω42μ2ω2q23q42(μ2ω23q2), (54)

and

Ω2=ω(μ2ω2+q2)ω(μ2ω2+q2)(μ2ω2+5q2)2(μ2ω2+q2). (55)

8  Stability Analysis for the Nonlinear Case

The first-order solvability condition (40) can be used to find the stability picture at the resonance case. The second-order solvability condition (42) can be used to find the value of the detuning parameter σ .

To relax the periodic term into the Eq. (40) we let

A(τ)=[α(τ)+iβ(τ)]eiστ, (56)

with real functions α and β . Insert Eq. (56) into Eq. (40), separating real and imaginary parts yields:

α+12μα(σ+q2ω+3k2ω(α2+β2))β=0, (57)

β+12μβ+(σq2ω3k2ω(α2+β2))α=0. (58)

In order to solve the above coupled nonlinear Eqs. (57) and (58), we may discuss the behavior at the steady-state response. This case is corresponding to the case of d..dτ=0 . If the solutions of Eqs. (57) and (58), at the steady-state, are represented by α0 and β0 , which are given by

12μα0(σ+q2ω+3kr22ω)β0=0, (59)

12μβ0+(σq2ω3kr22ω)α0=0, (60)

where r2=α02+β02 is used. Eqs. (59) and (60) are two coupled algebraic equations in α0 and β0 . For nontrivial solutions in α0 and β0 , we obtain

σ2=(q+3kr2)4ω2214μ2. (61)

Besides, the constants α0 and β0 may be chosen as

α0=(σ+q+3kr22ω), and β0=12μ. (62)

Squaring both equations in (62) and adding we get

(σ+q+3kr22ω)2=r214μ2. (63)

Combing Eq. (61) with Eq. (63) yields

σ=2ω2r2(q+3kr2)22ω(q+3kr2). (64)

In order to find a constrain for a bounded solution we may modulate the functions α and β as

α(τ)=α0+α1(τ)andβ(τ)=β0+β1(τ), (65)

where the functions α1(τ) and β1(τ) refer to a small deviation from the steady-state solution α0 and β0 . Then the system of Eqs. (57) and (58) in the linearizing form becomes

α1+(12μ3kωα0β0)α1(σ+q+3kr22ω+3kωβ02)β1=0, (66)

β+(12μ3kωα0β0)β1+(σq+3kr22ω3kωα02)α1=0. (67)

The above system is two coupled linear differential equations of first-order in the two functions α1 and β1 . This system can be satisfied by

α1(τ)=(σ+q+3kr22ω+3kωβ02)e(12μ3kωα0β0)τsinΘτ, (68)

β1(τ)=e(12μ3kωα0β0)τΘcosΘτ, (69)

where Θ is given by the following characteristic equation:

Θ2=(σq+9kr22ω+3k4ωμ2)(σ+q+3kr22ω+3k4ωμ2). (70)

where relations (62) are used. This characteristic equation depends on the two related parameters σ and r2 . This relation between them is given in Eq. (61) or in Eq. (64).

With the help of the second-order solvability condition (42) one can find an expression for both the unknowns σ and r2 in terms of the frequency Ω . To accomplish this, one may substitute the steady-state solution

A(τ)=(α0+iβ0)eiστ, (71)

into the second-order solvability condition (42). Separating the real and imaginary parts, produces the following relations, between the parameters σ , Ω and r2 :

σ2q24Ω(ω+Ω)+116kqμ2(Ω2+ωΩ6ω2ω2Ω(ω+Ω))[kq8(Ω2+ωΩ6ω2ω2Ω(ω+Ω))3k28ω2(σ+q+3kr22ω)]r2=0, (72)

σ+kq8(Ω2+ωΩ+6ω2ω2Ω(ω+Ω))(σ+q+3kr22ω)+3k22ω2r2=0, (73)

where relations (62) are used. Removing the parameter σ from Eq. (73), by using its equivalent in Eq. (64), gives a polynomial of second-order in r2 :

r4+Ω(ω+Ω)(12k2q+8ω323kqω)+6kqω336k2(kω)Ω(ω+Ω)r2q2ω9k2(kω)=0;kω. (74)

Replacing r4 and r2 into Eq. (72) with their equivalents in Eqs. (64) and (73) leads to the following quadratic equation in the detuning parameter σ :

6kωΩ(ω+Ω)(8ω21)[8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2)]σ2[Ω2(ω+Ω)2[2ω(8ω+kq)(3k2q3kq+4ω22ωq)+3k2q2]+6ω2Ω(ω+Ω)[3k2q2+8ω4+4kq2ω2+2kqω(3k2q3kq+4ω22ωq)]+144kq2ω6]σ+Ω2(ω+Ω)2[3k2q(ωμ2q)(8ω+kq)+kq2(3kq4ω2+2ωq)]36kω5(3k2μ2+2kq+2q3)+6kω2Ω(ω+Ω)[3k2q2(ωμ2q)qω(3kμ2+2q)(8ω+kq)+q2(3kq4ω2+2ωq)2q2ω]=0. (75)

This equation gives two values σ1 and σ2 for the detuning parameter σ which makes the solution (56) without unknowns.

The stabilization for the problem requires that the right-hand side of Eq. (70) be positive provided that the exponential in Eqs. (68) and (69) has positive values. It is noted that the stability reveals as the coefficient of the periodic term in Eq. (2) tends to zero. The instability arrived as the parameter q going away the zero value. Thus, the stability conditions are found as

μ>0,3kωσ+3kq+9k2r22ω21<0, (76)

σ>q+6kr22ω3k4ωμ2, (77)

σ+q+3kr22ω+3k4ωμ2<0. (78)

Removing the parameter r2 from the above stability conditions by using Eq. (73) yields the following two conditions for stability:

μ>0,ω(k1)σ<16(2ω23kq)+2kq2(Ω2+ωΩ+6ω2)8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2), (79)

3kσ>qω3q2(Ω2+ωΩ+6ω2)2ω[8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2)]. (80)

The transition curves separating the stable state of unstable one are corresponding to

σ=(2ω23kq)[8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2)]+12kq2(Ω2+ωΩ+6ω2)6ω(k1)[8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2)], (81)

and

σ=2kq[8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2)]3kq2(Ω2+ωΩ+6ω2)6ω[8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2)]. (82)

Using the definition (38) the above transition curves can be sought within the parameter ρ as

Ω=ω+ρ(2ω23kq)[8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2)]+12kq2(Ω2+ωΩ+6ω2)6ω(k1)[8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2)], (83)

Ω=ω+ρ2kq[8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2)]3kq2(Ω2+ωΩ+6ω2)6ω[8ωΩ(ω+Ω)+kq(Ω2+ωΩ+6ω2)]. (84)

To obtain the transition curves, independent of the parameter ρ , we may be inserting Eq. (81) as well as Eq. (82), into the relation (75), then the following transition curves are imposed

a3Ω3(Ω+ω)3+a2ω2Ω2(Ω+ω)2+36a1ω4Ω(Ω+ω)+a0=0, (85)

b3Ω3(Ω+ω)3+b2ω2Ω2(Ω+ω)2+36b1ω4Ω(Ω+ω)+b0=0. (86)

It is noted that the instability state lies between the above transition curves. The constant coefficients aj and bj, j=0,1,2,3 are given below:

a3=k(8ω21)[(2ω23kq)(8ω+kq)+12kq2]2(k1)[(2ω23kq)(8ω+kq)+12kq2][2ω(8ω+kq)(3k2q3kq+4ω22ωq)+3k2q2]+6ω(k1)2(8ω+kq)[3k2q(ωμ2q)(8ω+kq)+kq2(3kq4ω2+2ωq)],

a2=12k2q(8ω21)(2ω23kq+12q)[(2ω23kq)(8ω+kq)+12kq2]6(k1)[(2ω23kq)(8ω+kq)+12kq2][3k2q2+8ω4+4kq2ω2+2kqω(3k2q3kq+4ω22ωq)]6kq(k1)(2ω23kq+12q)[2ω(8ω+kq)(3k2q3kq+4ω22ωq)+3k2q2]+6×36kω(k1)2k2q2(ωμ2q)(8ω+kq)+36kω(k1)2q2(3kq4ω2+2ωq)(8ω+2kq)36kω(k1)2qω(3kμ2+2q)(8ω+kq)272kω(k1)2q2ω(8ω+kq),

a1=+k3q2(8ω21)(2ω23kq+12q)24kq2ω2(k1)[(2ω23kq)(8ω+kq)+12kq2]kq(k1)(2ω23kq+12q)[3k2q2+8ω4+4kq2ω2+2kqω(3k2q3kq+4ω22ωq)]+6k2qω(k1)2[3k2q2(ωμ2q)qω(3kμ2+2q)(8ω+kq)+q2(3kq4ω2+2ωq)2q2ω]6kω2(k1)2(8ω+kq)(3k2μ2+2kq+2q3),

a0=144kq2ω56kqω3(k1)(2ω23kq+12q)36×36k2qω8(k1)2(3k2μ2+2kq+2q3),

b3=k2q2(8ω21)[2(8ω+kq)3q]2+6qω(8ω+kq)[3k(ωμ2q)(8ω+kq)+q(3kq4ω2+2ωq)][2q(8ω+kq)3q2][2ω(8ω+kq)(3k2q3kq+4ω22ωq)+3k2q2],

b2=12k2q3(8ω21)(2k3)[2(8ω+kq)3q]+6kq2[3k(ωμ2q)(8ω+kq)+q(3kq4ω2+2ωq)]6q[2(8ω+kq)3q][3k2q2+8ω4+4kq2ω2+2kqω(3k2q3kq+4ω22ωq)]6q2(2k3)[2ω(8ω+kq)(3k2q3kq+4ω22ωq)+3k2q2]+36ω(8ω+kq)[3k2q2(ωμ2q)qω(3kμ2+2q)(8ω+kq)+q2(3kq4ω2+2ωq)2q2ω],

b1=k2q4(8ω21)(2k3)24q2ω2[2kq(8ω+kq)3kq2]6ω2(3k2μ2+2kq+2q3)(8ω+kq)q(2k3)[3k2q2+8ω4+4kq2ω2+2kqω(3k2q3kq+4ω22ωq)]+6kqω[3k2q2(ωμ2q)qω(3kμ2+2q)(8ω+kq)+q2(3kq4ω2+2ωq)2q2ω],

b0=72×6ω8kq(9k2μ2+6kq+4kq3).

9  Conclusion

The homotopy perturbation method (HPM) is one of a easy, powerful, efficient, and accurate approach for evalueting solutions of a large class of nonlinear equations without the need of a discretization or linearization process. HPM is a combination of the homotopy and perturbation methods. That can take the advantages of the conventional perturbation method and eliminating its restrictions. It yields a rapid convergence of the solution series with a few iterations leading to accurate solutions, and the round-off errors are avoided. In general, this method has been successfully used to solve different kinds of linear and nonlinear problems in engineering and science. So, in the present work, we propose a variation of the homotopy perturbation approach via a modulation method that allows finding analytic solutions for ordinary differential models with periodic coefficients. This article is prepared to analyze a parametrically excited oscillator in the presence of strong cubic nonlinearity. The simplest model of this kind is the Mathieu equation that contains a small parameter [47,48]. The present analysis that employs the homotopy perturbation approach , has no dependence on models having a small parameter. Due to the present modulation approach, at each level of perturbation, a solvability condition is enjoined. By solving these solvability conditions drives to examining the stability behavior. In each resonance/non-resonant cases stability conditions, are obtained.

Author Contributions: The first author proposed and developed the mathematical modeling of the problem and examines the theory validation; the second author analyzed the cubic damping Mathieu equation. The manuscript was written through the contribution of all authors. All authors discussed the outcomes, reviewed, and approved the final version of the manuscript.

Funding Statement: The authors received no financial support for this research, authorship and publication of this article.

Conflicts of Interest: The authors declare that there are no competing interests regarding the publication of the present paper.

## References

1. Rong, H., Xu, W., & Fang, T. (1998). Principle response of Duffing oscillator to combined deterministic and narrow-band random parametric excitation. Journal of Sound and Vibration, 210(4), 483-515. [Google Scholar] [CrossRef]
2. Morrison, T. M. (2016). Three problems in nonlinear dynamics with 2:1 parametric excitation (Ph.D. thesis). Department of Applied Mechanics, University of Glasgow, Scotland, UK.
3. Dutta, T. K., & Prajapati, P. K. (2016). Some dynamical properties of the Duffing equation. International Journal of Engineering Research and Technology, 5(12), 500-503. [Google Scholar] [CrossRef]
4. Al-Jawary, M. A., & Abd-Al-Razaq, S. G. (2016). Analytic and numerical solution for Duffing equations. International Journal of Basic and Applied Sciences, 5(2), 115-119. [Google Scholar] [CrossRef]
5. Sunday, J. (2017). The Duffing oscillator: Applications and computational simulations. Asian Research Journal of Mathematics, 2(3), 1-13. [Google Scholar] [CrossRef]
6. Mclachlan, N. W. (1951). Theory and application of Mathieu functions. London, UK: Oxford University Press.
7. Ramani, D. V., Keith, W. L., & Rand, R. H. (2004). Perturbation solution for secondary bifurcation in the quadratically-damped Mathieu equation. International Journal of Non-Linear Mechanics, 39(3), 491-502. [Google Scholar] [CrossRef]
8. Taylor, J. H., & Narendra, K. S. (1969). Stability regions for the damped Mathieu equation. SIAM Journal on Applied Mathematics, 17(2), 343-352. [Google Scholar] [CrossRef]
9. Insperger, T., & Stépán, G. (2003). Stability of the damped Mathieu equation with time delay. Journal of Dynamic Systems, Measurement, and Control, 125(2), 166-171. [Google Scholar] [CrossRef]
10. El-Nady, A. O., & Lashin, M. M. A. (2016). Approximate solution of nonlinear Duffing oscillator using Taylor expansion. Journal of Mechanical Engineering and Automation, 6(5), 110-116. [Google Scholar] [CrossRef]
11. Feng, Z., Chen, G., & Hsu, S. B. (2006). A qualitative study of the damped Duffing equation and applications. Discrete and Continuous Dynamical Systems–Series B, 5(5), 1-24. [Google Scholar] [CrossRef]
12. Luo, A. C. J., & Ma, H. (2018). Bifurcation trees of periodic motions to chaos in a parametric Duffing oscillator. International Journal of Dynamics and Control, 6(2), 425-458. [Google Scholar] [CrossRef]
13. Michon, G., Manin, L., Parker, R. G., & Dufour, R. (2008). Duffing oscillator with parametric excitation: Analytical and experimental investigation on a belt-pulley system. Journal of Computational and Nonlinear Dynamics, 3(3), 31001. [Google Scholar] [CrossRef]
14. Zivieri, R., Vergura, S., & Carpentier, M. (2016). Analytical and numerical solution to the nonlinear cubic Duffng equation: An application to electrical signal analysis of distribution lines. Applied Mathematical Modelling, 40(21–22), 9152-9164. [Google Scholar] [CrossRef]
15. Moatimid, G. M. (2020). Stability analysis of a parametric duffing oscillator. Journal of Engineering Mechanics, 146(5), 5020001. [Google Scholar] [CrossRef]
16. Nayfeh, A. H., Mook, D. T. (1979). Non-linear oscillations. New York, USA: John Wily & Sons.
17. He, J. H. (1999). Homotopy perturbation technique. Computational Methods in Applied Mechanics and Engineering, 178(3–4), 257-262. [Google Scholar] [CrossRef]
18. He, J. H. (2000). A coupling method of homotopy technique and a perturbation technique for non-linear problems. International Journal of Non-Linear Mechanics, 35(1), 37-43. [Google Scholar] [CrossRef]
19. He, J. H. (2005). Application of homotopy perturbation method to nonlinear wave equations. Chaos, Solitons and Fractals, 26(3), 295-300. [Google Scholar] [CrossRef]
20. He, J. H. (2006). Homotopy perturbation method for solving boundary value problems. Physics Letters A, 350(1–2), 87-88. [Google Scholar] [CrossRef]
21. He, J. H. (2012). Homotopy perturbation method with an auxiliary term. Abstract and Applied Analysis, 2012, [Google Scholar] [CrossRef]
22. He, J. H. (2003). Homotopy perturbation method: A new nonlinear analytical technique. Applied Mathematics and Computation, 135(1), 73-79. [Google Scholar] [CrossRef]
23. He, J. H. (2003). A simple perturbation approach to Blasius equation. Applied Mathematics and Computation, 140(2–3), 217-222. [Google Scholar] [CrossRef]
24. Babolian, E., Saeidian, J., & Azizi, A. (2009). Application of homotopy perturbation method to some nonlinear problems. Applied Mathematical Sciences, 3(45), 2215-2226. [Google Scholar]
25. El-Dib, Y. O., & Moatimid, G. M. (2019). Stability configuration of a rocking rigid rod over a circular surface using the homotopy perturbation method and Laplace transform. Arabian Journal for Science and Engineering, 44(7), 6581-6659. [Google Scholar] [CrossRef]
26. He, J. H. (2014). Homotopy perturbation method with two expanding parameters. Indian Journal of Physics, 88(2), 193-196. [Google Scholar] [CrossRef]
27. El-Dib, Y. O. (2018). Multi-homotopy perturbations technique for solving nonlinear partial differential equations with Laplace transforms. Nonlinear Science Letters A, 9, 349-359. [Google Scholar]
28. El-Dib, Y. O. (2017). Multiple scales homotopy perturbation method for nonlinear oscillators. Nonlinear Science Letters A, 8(4), 352-364. [Google Scholar]
29. Ren, Z. F., Yao, S. W., & He, J. H. (2019). He’s multiple scales method for nonlinear vibrations. Journal of Low Frequency Noise, Vibration & Active Control, 38(3–4), 1708-1712. [Google Scholar] [CrossRef]
30. Rabbani, M. (2013). New homotopy perturbation method to solve non-linear problems. Journal of Mathematics and Computer Science, 7(4), 272-275. [Google Scholar] [CrossRef]
31. Ayati, Z., Biazar, J., & Ebrahimi, S. (2014). A new homotopy perturbation method for solving linear and nonlinear Schrödinger equations. Journal of Interpolation and Approximation in Scientific Computing, 2014(1), 1-8. [Google Scholar] [CrossRef]
32. El-Dib, Y. O., Moatimid, G. M., & Mady, A. A. (2019). A novelty to the nonlinear rotating Rayleigh-Taylor instability. Pramana–Journal of Physics, 93, 82. [Google Scholar] [CrossRef]
33. El-Dib, Y. O. (2019). Periodic solution of the cubic nonlinear Klein-Gordon equation and the stability criteria via the He-multiple-scales method. Pramana–Journal of Physics, 92, 7. [Google Scholar] [CrossRef]
34. El-Dib, Y. O., & Elgazery, N. S. (2020). Effect of fractional derivative properties on the periodic solution of the nonlinear oscillations, Fractals. Fractals, 28(7), 2050095. [Google Scholar] [CrossRef]
35. El-Dib, Y. O., & Mady, A. A. (2020). He’s multiple-scale solution for the three-dimensional nonlinear KH instability of rotating magnetic fluids. International Annals of Science, 9(1), 52-69. [Google Scholar] [CrossRef]
36. El-Dib, Y. O. (2020). Stability approach of a fractional-delayed duffing oscillator. Discontinuity Nonlinearity and Complexity, 9(3), 367-376. [Google Scholar] [CrossRef]
37. Moatimid, G. M., El-Dib, Y. O., & Zekry, M. H. (2020). The nonlinear instability of a cylindrical interface between two hydromagnetic darcian flows. Arabian Journal for Science and Engineering, 45(1), 391-409. [Google Scholar] [CrossRef]
38. El-Dib, Y. O. (2020). Modified multiple scale technique for the stability of the fractional delayed nonlinear oscillator. Pramana–Journal of Physics, 94(1), 56. [Google Scholar] [CrossRef]
39. Shen, Y., & El-Dib, Y. O. (2020). A periodic solution of the fractional sine-Gordon equation arising in architectural engineering. Journal of Low Frequency Noise, Vibration and Active Control, 40(2), 683-691. [Google Scholar] [CrossRef]
40. Elgazery, N. S. (2020). A periodic solution of the newell-whitehead-segel (NWS) wave equation via fractional calculus. Journal of Applied and Computational Mechanics, 6, 1293-1300. [Google Scholar] [CrossRef]
41. Yu, D. N., He, J. H., & Garcıa, A. G. (2019). Homotopy perturbation method with an auxiliary parameter for nonlinear oscillators. Journal of Low Frequency Noise, Vibration and Active Control, 38(3–4), 1540-1554. [Google Scholar] [CrossRef]
42. He, J. H., & El-Dib, Y. O. (2020). Homotopy perturbation method for Fangzhu oscillator. Journal of Mathematical Chemistry, 58(10), 2245-2253. [Google Scholar] [CrossRef]
43. He, J. H., & Jin, X. (2020). A short review on analytical methods for the capillary oscillator in a nanoscale deformable tube. Mathematical Methods in the Applied Sciences, 38(3), 1676. [Google Scholar] [CrossRef]
44. He, C. H., He, J. H., Sedighi, H. M. (2020). Fangzhu (方诸): An ancient Chinese nanotechnology for water collection from air: History, mathematical insight, promises and challenges. Mathematical Methods in the Applied Sciences,2020, 1–10. DOI 10.1002/mma.6384. [CrossRef]
45. He, J. H. (2020). A short review on analytical methods for to a fully fourth-order nonlinear integral boundary value problem with fractal derivatives. International Journal of Numerical Methods for Heat and Fluid Flow, 30(11), 4933-4943. [Google Scholar] [CrossRef]
46. Nash, C., Sen, S. (1983). Topology and geometry for physicists. London, UK: Academic Press.
47. El-Dib, Y. O. (2001). Nonlinear Mathieu equation and coupled resonance mechanism. Chaos, Solitons & Fractals, 12(4), 705-720. [Google Scholar] [CrossRef]
48. El-Dib, Y. O. (2016). Stability criterion for time-delay 3-dimension damped Mathieu equation. Science and Engineering Applications, 1(5), 76-88. [Google Scholar] This work is licensed under a Creative Commons Attribution 4.0 International License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.