Abstract

The objective of this paper is to present an analytical investigation to analyze the vibration of parametrically excited oscillator with strong cubic negative nonlinearity based on Mathieu-Duffing equation. The analytic investigation was conducted by using He's homotopy-perturbation method (HPM). In order to obtain the analytical solution of Mathieu-Duffing equation, homotopy-perturbation method has been utilized. The Runge-Kutta's (RK) algorithm was used to solve the governing equation via numerical solution. Finally, to demonstrate the validity of the proposed method, the response of the oscillator, which was obtained from approximate solution, has been shown graphically and compared with that of numerical solution. Afterward, the effects of variation of the parameters on the accuracy of the homotopy- perturbation method were studied.

1. Introduction

Parametrically excited systems are widely spread in many branches of physics and engineering. It is important to investigate the dynamic behavior of these systems.

Vibrations of such systems occur in wide range of mechanics, due to time-varying loads, especially periodic ones. These vibrations appear in columns made of nonlinear elastic material [1], beams with a harmonically variable length [2], beams with harmonic motion of their support [3], floating offshore structures [4], parametrically excited pendulums [5], cables being towed by a submarine [6, 7], and so forth. Parametric excitations occur in electrostatically driven microelectromechanical oscillators [8], which are produced by fluctuating voltages applied across comb drives. In practical engineering situations, the properties of parametric oscillations are widely used, for example, in the radio, the computer, and laser engineering, in vibrant machines with special design [9], Paul trap mass spectrometers [10], as well as in a simulator for proving the equivalence of inertia and passive gravitational mass [11]. Parametric resonance has been well established in many areas of science, including the stability of ships, the forced motion of a swing, and Faraday surface wave patterns on water. The highly sensitive mass sensor is studied as an in-plane parametrically resonant oscillator [12]. The simplest mathematical model of the system with a parametric periodic load is usually a linear Mathieu differential equation. Due to the nonlinear properties of a real system, nonlinear terms are added to the equation by Younesian et al. [13]. They have added a cubic nonlinearity term (π‘₯3) to Mathieu-Duffing equation using averaging method [13, 14]. It is obvious that for such nonlinear equations, there is no precise analytical solution; consequently these nonlinear equations need to be solved using other methods. In recent decades, numerical methods have been well established for analyzing the nonlinear equations such as Mathieu-differential equation. Most scientists believe that the combination of numerical and semianalytical methods can give rise to useful results. The two common methods employed to solve the Mathieu-differential nonlinear equation are based on the numerical integration by Jazar [15] and perturbation by Cveticanin and Kovacic [16]. Perturbation method, which is much more developed, is one of the well-known methods that were studied by some researchers [17–19] to solve the linear and nonlinear equations. Actually, these scientists had paid more attention to the mathematical aspects. Since there are some limitations in using the common perturbation method together with the fact that this method is based on the existence of small parameter, developing the method for different applications is very difficult. Therefore, a number of semianalytical methods, including the variational iteration method [20, 21] and the homotopy-perturbation method (HPM), [22–24] have recently been developed by using new techniques to eliminate the small parameter. Other authors have made much contribution to the applications of these methods [25–37].

In this paper, we extend the work of Cveticanin and Kovacic [16] by applying HPM to analyze the vibrations of parametrically excited oscillator with strong cubic negative nonlinearity. At the first, a brief review of the method will be given and then the HPM will be applied to the Mathieu-Duffing equation (1.1). Next, Runge-Kutta’s (RK) algorithm will be introduced to solve the governing equation. Finally, a numerical example is given to demonstrate the validity of the proposed method (HPM) and the effect of the variation parameters on the accuracy of the HPM is investigated.

The governing equation of Mathieu-Duffing system which is considered in this study is described by the following high-order nonlinear differential equation:

ξ€Ίξ€·Μˆπ‘₯+𝛿+2πœ€cos2𝑑π‘₯βˆ’πœ™π‘₯3=0,(1.1)where dots indicate differentiation with respect to the time (t), πœ€β‰ͺ1 is a small parameter, πœ‘ is the parameter of nonlinearity, and Ξ΄ is the transient curve and can be defined as [16]

𝛿=πœ™π‘₯20ξ‚€1βˆ’2πœ€2+πœ™π‘₯20.(1.2)The initial condition considered in this study is defined by π‘₯ = 0.1, Μ‡π‘₯ = 0 [16].

This mathematical model corresponds to the parametrically excited oscillator with a softening spring.

2. Analytical Solution

The analytical solution of such high-order nonlinear differential equation as (1.1) is usually difficult to find. For these equations, much attention was paid to develop various simple numerical and analytical methods such as those based on the numerical integration [15] and perturbation [16]. In this paper, an approximate analytical method is developed to solve the problem. However, as mentioned before, for this case, some analytical solutions are well established. Application of the proposed semianalytical method for achieving the nonlinearity of parametrically excited vibrations of an oscillator is also discussed in this paper.

2.1. Homotopy-Perturbation Method (HPM)
2.1.1. Background

In this paper, we apply HPM [22–37] for the solution of the Mathieu-Duffing equation (1.1). In order to demonstrate how this method works, let us consider the following nonlinear differential equation:

𝐴(𝑒)βˆ’π‘“(π‘Ÿ)=0,π‘ŸβˆˆΞ©,(2.1)subject to the boundary conditions of

𝐡𝑒,πœ•π‘’ξ‚πœ•π‘Ÿ=0,π‘ŸβˆˆΞ“,(2.2)where A is a general differential operator, B is a boundary operator, 𝑓(π‘Ÿ) is a known analytic functional, and Ξ“ is the boundary of the domain Ξ©.

The operator A can generally be divided into two parts L and N, where L is linear, whereas N is nonlinear. Therefore, (2.1) can be rewritten as

𝐿(𝑒)+𝑁(𝑒)βˆ’π‘“(π‘Ÿ)=0.(2.3)Homotopy-perturbation structure is shown to be as the following equation [22–37]:

𝐻(𝑣,𝑝)=(1βˆ’π‘)𝐿(𝑣)βˆ’πΏ(𝑒0)ξ€»+𝑝[𝐿(𝑣)+𝑁(𝑣)βˆ’π‘“(π‘Ÿ)]=0,(2.4)where

𝑣(π‘Ÿ,𝑝)βˆΆΞ©Γ—[0,1]→𝑅.(2.5)For p = 0 and p = 1, (2.4) reduces to the following equations, respectively:

𝐻(𝑣,0)=𝐿(𝑣)βˆ’πΏ(𝑒0)=0,𝐻(𝑣,1)=𝐴(𝑣)βˆ’π‘“(π‘Ÿ)=0,(2.6)where π‘βˆˆ[0,1] is an embedding parameter and 𝑒0 is the first approximation that satisfied the boundary condition.

The process of changes in p from zero to unity is that of 𝑣(π‘Ÿ,𝑝) changing form 𝑒0 to 𝑒(π‘Ÿ). We consider v, as the following [22–37]:

𝑣=𝑣0+𝑝𝑣1+𝑝2𝑣2+β‹―,(2.7)and the best approximation for solution is [22–37]𝑒=lim𝑝→1𝑣=𝑣0+𝑣1+𝑣2+β‹―.(2.8)

2.1.2. Implementation of HPM

Now we apply homotopy-perturbation (2.4) to (1.1). We construct a homotopy in the following form:

𝑑(1βˆ’π‘)2π‘₯(𝑑)𝑑𝑑2ξ‚„ξ€½βˆ’ξ€Ίξ€·π‘₯𝑑𝑑+𝑝𝛿+2πœ€cos2𝑑+πœ™π‘₯3ξ€Ύ=0.(2.9)Obviously, (2.9) becomes linear if p = 0, when p = 1, it becomes the original nonlinear one. So the variation of p from zero to unity makes the equation to change to nonlinear from linear one.

According to HPM, we assume that the solution of (2.9) can be expressed in a series of p

𝑣(π‘Ÿ)=𝑣0(π‘Ÿ)+𝑝𝑣1(π‘Ÿ)+𝑝2𝑣2(π‘Ÿ)+β‹―.(2.10)Substituting 𝑣(π‘Ÿ) from (2.10) into (2.9), after some simplification and substitution and rearranging based on powers of p-terms, we have

𝑝0βˆΆπ‘‘2𝑣(𝑑)𝑑𝑑2𝑝=0,(2.11)1βˆΆπ‘‘2𝑣1(𝑑)𝑑𝑑2+𝛿+4πœ€cos(𝑑)2βˆ’πœ™π‘£0(𝑑)2ξ€»π‘£βˆ’2πœ€0𝑝(𝑑)=0,(2.12)2βˆΆπœ€4Γ—103ξ€Ίβˆ’102ξ€·ξ€Έπ‘‘πœ€+πœ‘βˆ’100𝛿2+2Γ—102ξ€»βˆ’βˆ’3πœ‘cos(2𝑑)π‘‘πœ€2Γ—103ξ€·πœ‘βˆ’102π›Ώξ€Έπœ€sin(2𝑑)+2cos(4𝑑)3.2Γ—102+ξ€·2.1875Γ—10βˆ’2βˆ’2.5Γ—10βˆ’2𝑑2ξ€Έπœ€2+πœ€24Γ—105ξ€Ίξ€·βˆ’1.8Γ—103πœ‘+6Γ—104𝛿𝑑2+1.8Γ—103πœ‘βˆ’1.2Γ—105𝛿+ξ€·πœ‘βˆ’102𝛿/3ξ€Έξ€·πœ‘βˆ’102𝛿𝑑48Γ—105.(2.13)Here, the initial condition at this boundary can be determined by the boundary condition (see Section 1).

With the effective initial approximation for 𝑒0 from the initial conditions to (2.11), we construct 𝑒0:

𝑣0(𝑑)=0.1.(2.14)Then, solving (2.12), (2.13), we have

𝑣1𝑑=ξ‚€πœ™βˆ’π›Ώ2002𝑑2βˆ’πœ€2+πœ€cos(2𝑑)2,𝑣21(𝑑)=4Γ—103πœ€ξ€Ίβˆ’102ξ€·πœ€+βˆ’102𝑑𝛿+πœ‘2+2Γ—102ξ€»ξ€·ξ€Έβˆ’1π›Ώβˆ’3πœ‘cos2𝑑2Γ—103ξ€·πœ€π‘‘βˆ’102ξ€Έξ€·ξ€Έ+1𝛿+πœ‘sin2𝑑3.2Γ—102πœ€2ξ€·ξ€Έ+1cos4𝑑2.4Γ—106ξ€·5.25Γ—104βˆ’6Γ—104𝑑2ξ€Έπœ€2+12.4Γ—106ξ€Ίξ€·βˆ’1.8Γ—103πœ‘+6Γ—104𝛿𝑑2βˆ’1.2Γ—105𝛿+1.8Γ—103πœ‘ξ€»πœ€+18Γ—105ξ€·βˆ’102𝛿+πœ‘πœ‘βˆ’1023𝛿𝑑4.(2.15)So substituting (2.14)-(2.15) into (2.10) gives the solution in the following form:

πœ‘π‘£(𝑑)=0.1+𝑝2Γ—102βˆ’π›Ώ2𝑑2βˆ’πœ€2+πœ€cos(2𝑑)2ξ‚„+𝑝214Γ—103πœ€ξ€Ίβˆ’102ξ€·πœ€+βˆ’102𝑑𝛿+πœ‘2+2Γ—102ξ€»ξ€·ξ€Έβˆ’1π›Ώβˆ’3πœ“cos2𝑑2Γ—103ξ€·πœ€π‘‘βˆ’102ξ€Έξ€·ξ€Έ+1𝛿+πœ‘sin2𝑑3.2Γ—102πœ€2ξ€·ξ€Έ+1cos4𝑑2.4Γ—106ξ€·5.25Γ—104βˆ’6Γ—104𝑑2ξ€Έπœ€2+12.4Γ—106ξ€Ίξ€·βˆ’1.8Γ—103πœ‘+6Γ—104𝛿𝑑2βˆ’1.2Γ—105𝛿+1.8Γ—103πœ‘ξ€»πœ€+18Γ—105ξ€·βˆ’102𝛿+πœ‘πœ‘βˆ’1023𝛿𝑑4+β‹―.(2.16)As 𝑝→1,𝑣(𝑑)β†’π‘₯(𝑑) in (2.16), the solution for the displacement function of vibrations in a parametrically excited oscillator with strong cubic negative nonlinearity is given. From (2.14)-(2.15), it can be realized that HPM requires only three steps to get more accurate results rather than other analytical methods.

3. Numerical Solution

In order to provide a basis for the purpose of comparison of analysis of parametrically excited vibrations of an oscillator with strong cubic negative nonlinearity with no explicit solution, in this section, the solution of the governing equation of the system (1.1) by using RK iterative method is presented. The method can also be used for other second-order differential equations.

3.1. Runge-Kutta's Method (RK)

As shown by (1.1), the solution of vibration of parametrically excited oscillator with strong cubic negative nonlinearity is to solve a second-order differential equation with a set of certain boundary conditions. For this case, (1.1) can be written in the following general form:

ξ€·ξ€ΈΜˆπ‘₯=𝑓(𝑑,π‘₯,Μ‡π‘₯)=βˆ’π›Ώ+2πœ€cos2𝑑π‘₯+πœ™π‘₯.(3.1)For such a boundary value problem given by boundary condition defined in Section 4, some numerical methods have been developed [38, 39]. Here we apply the fourth-order RK algorithm to solve (3.1) subject to the given boundary conditions. RK iterative formulae for the second-order differential equation are [39]

Μ‡π‘₯𝑖+1=Μ‡π‘₯𝑖+Δ𝑑6ξ€·π‘˜1+2π‘˜2+2π‘˜3+π‘˜4ξ€Έ,π‘₯𝑖+1=π‘₯𝑖+Δ𝑑̇π‘₯𝑖+Δ𝑑6ξ€·π‘˜1+π‘˜2+π‘˜3ξ€Έξ‚„,(3.2)where Δ𝑑 is the increment of the time and π‘˜1,π‘˜2,π‘˜3, and π‘˜4 are determined from the following formulae:

π‘˜1𝑑=𝑓𝑖,π‘₯𝑖,Μ‡π‘₯𝑖,π‘˜2𝑑=𝑓𝑖+Δ𝑑2,π‘₯𝑖+Δ𝑑2Μ‡π‘₯𝑖,Μ‡π‘₯𝑖+Δ𝑑2π‘˜1,π‘˜3𝑑=𝑓𝑖+Δ𝑑2,π‘₯𝑖+Δ𝑑2Μ‡π‘₯𝑖+14Δ𝑑2π‘˜1,Μ‡π‘₯𝑖+Δ𝑑2π‘˜2,π‘˜4𝑑=𝑓𝑖+Δ𝑑,π‘₯𝑖+Δ𝑑̇π‘₯𝑖+12Δ𝑑2π‘˜2,Μ‡π‘₯𝑖+Ξ”π‘‘π‘˜3.(3.3)The numerical solution starts from the boundary at the initial time, where the first value of the displacement function and its first-order derivative is determined from initial condition (see Section 1). Then, with a small time increment (Δ𝑑), the displacement function and its first-order derivative at the new position can be obtained using (3.2). This process continues to the end of time.

4. Results and Discussion

In this study, the Mathieu-Duffing equation has been solved by using HPM and RK method. The results shown in Figures 1–6 indicate that the HPM experiences a high accuracy. In addition, in comparison with RK method and other analytical methods, a considerable reduction of the volume of the calculation can be seen in HPM. It can be approved that HPM is powerful in finding analytical solutions for a wide class of nonlinear problems.

Figures 1(a) and 1(b) illustrate the time history diagram of the displacement and velocity, respectively. Figure 2 illustrates the velocity versus the displacement obtained from the proposed method. Figure 1 as well as Figure 2 are obtained for πœ‘ = 2, Ξ΅ = 0.01. These figures show obviously the excellent agreement between HPM and RK method. Figure 2 also shows the residual loop; the area of loops corresponds to the loss of energy in each cycle. It illustrates that the increase of cycle numbers reduces the loss of energy.

For further verification, the accuracy of Figure 1 is shown in Figure 3 which shows that increasing the time gives rise to reduce the accuracy of displacement as well as the velocity. Figure 4 confirms the accuracy of results depicted in Figure 2. Figure 5 shows the effects of the variation of the parameter Ξ΅ at t = 7 and πœ‘ = 2 on the accuracy of HPM. From this figure, it is clear that the change of the Ξ΅ from 0.01 to 1 reduces the accuracy up to 10%. So, it can be seen that the best approximation for high accurate results can be obtained in the range of Ξ΅ less than 0.1. From Figure 6, it is evident that at t = 7 and Ξ΅ = 0.01, the variation of πœ‘ from 2 to 4 reduces the accuracy up to 95%.

5. Conclusion

In this survey, the HPM has been employed to analyze the parametrically excited vibrations of an oscillator with strong cubic negative nonlinearity. The results obtained from this method have been compared with those obtained from numerical method using RK algorithm. This comparison shows excellent agreement between the two methods.

The presented scheme provides concise and straightforward solution to approach reliable results, and it overcomes the difficulties that have been arisen in conventional methods. Also, HPM does not require small parameters, so the limitation of the conventional perturbation method could be eliminated.

Solution of the Mathieu-Duffing equation shows that accuracy of the results is affected by the variation of the parameters Ξ΅ and πœ‘ considerably. In other words, increasing of the Ξ΅ and πœ‘ results in a reduction in the accuracy.

Nomenclature
t:Time
Μ‡π‘₯:Displacement
Μ‡π‘₯:Velocity
Ξ΅:Small parameter
Ο†:Parameter of nonlinearity
Ξ΄:Transient curve
p:Homotopy-perturbation parameter