Abstract
The dynamical characteristics of Mathieu equation with fractionalorder derivative is analytically studied by the LindstedtPoincare method and the multiplescale method. The stability boundaries and the corresponding periodic solutions on these boundaries for the constant stiffness ${\delta}_{0}={n}^{2}$ ($n$ = 0, 1, 2, …), are analytically obtained. The effects of the fractionalorder parameters on the stability boundaries and the corresponding periodic solutions, including the fractional coefficient and the fractional order, are characterized by the equivalent linear damping coefficient (ELDC) and the equivalent linear stiffness coefficient (ELSC). The comparisons between the transition curves on the boundaries obtained by the approximate analytical solution and the numerical method verify the correctness and satisfactory precision of the analytical solution. The following analysis is focused on the effects of the fractional parameters on the stability boundaries located in the $\delta \epsilon $ plane. It is found that the increase of the fractional order $p$ could make the ELDC larger and ELSC smaller, which could result into the rightwards and upwards moving of the stability boundaries simultaneously. It could also be concluded the increase of the fractional coefficient ${K}_{1}$ would make the ELDC and ELSC larger, which could move the transition curves to the left and upwards at the same time. These results are very helpful to design, analyze or control this kind of system, and could present beneficial reference to the similar fractionalorder system.
1. Introduction
In 1695, the French mathematician Hospital and German mathematician Leibniz put forward the concept of fractionalorder calculus for the first time. In the following centuries, the theory about fractionalorder calculus including the definition, characteristics, calculation of fractionalorder calculus and the relationship with integerorder calculus, etc, had been developing rapidly and a series of meaningful results had been obtained [17].
In recent years, fractionalorder calculus was paid more and more attention from researchers in different fields and became an international hot topic. A lot of mathematicians, physicists, chemists, dynamicists, and engineers in some relevant fields had applied this mathematical tool to solve the problems they met. For example, Shen and Yang et al. [811] investigated several linear and nonlinear fractionalorder oscillators by averaging method, and found that the fractionalorder derivatives had both damping and stiffness effects on the dynamical response in those oscillators. Gorenflo et al. [12], Jumarie [13], Ishteva et al. [14] and Agnieszka et al. [15] respectively studied the definitions and numerical methods of fractionalorder calculus for GrünwaldLetnikov, RiemannLiouville and Caputo. Chen and Zhu [16] reviewed the analytical methods and control strategies for quasiintegrable Hamiltonian systems with fractionalorder derivatives, and pointed out some possible developing issues. They [1720] also studied some nonlinear fractionalorder system with different kinds of noise, and obtained some important statistic properties of the fractionalorder system. Wang et al. [2122] investigated a linear single degreeoffreedom oscillator with fractionalorder derivative, and obtained the composition solution of its initial value problem without external excitation. Li et al. [23] discussed the properties of three kinds of fractional derivatives and sequential property of the Caputo derivative is also derived. Li et al. [2426] had done a lot of researches in the mathematical theory of fractionalorder calculus, and also established some efficient numerical algorithms. Wahi and Chatterjee [27] studied a special linear single degreeoffreedom oscillator with fractionalorder derivative by averaging method, and analyzed the effects of the fractionalorders derivative. Xu et al. [28] used LindstedtPoincaré method and the multiplescale approach to investigate fractionalorder Duffing oscillator subjected to random excitation.
In addition to the theoretical research, fractional calculus had also been used to solve engineering problems. Comparing with the traditional integerorder system, the fractionalorder system is much closer to the real nature of the world, and has more advantages, such as strong ability of antinoise, good robustness, high control precision and so on. Accordingly fractionalorder calculus has the significant value for the engineering field.
The wellknown Mathieu equation, first introduced by Mathieu [29] when he studied the vibration of elliptical membranes in 1868, is a linear differential equation with periodic coefficients. Mathieu equation had been applied in physics and engineering fields, and many complex dynamical properties had been found herein. In recent years, many scholars had studied the Mathieu equation and found that the fractionalorder system could generate different dynamical properties from the integerorder system. For example, the chaotic behaviors of a nonlinear damped Mathieu system were investigated by applying numerical integration method [30]. Ebaid et al. [31] considered the fractionalorder calculus model of damped Mathieu equation and obtained the approximate analytical solution by using two different methods. The transition curves separating the stable and unstable regions in the fractionalorder Mathieu equation had been studied by the method of harmonic balance [3233], and Leung et al. also investigated a general version of fractionalorder MathieuDuffing equation by harmonic balance method.
In this paper, the fractionalorder Mathieu equation with linear damping is studied by the LindstedtPoincaré (LP) method coupled with the multiplescale method, and the effects of the fractional parameters on the dynamical properties of the fractionalorder Mathieu equation are analyzed in detail. The paper is organized as follow. Section 2 presents the stability boundaries and the corresponding approximately analytical periodic solutions on these boundaries for the constant stiffness ${\delta}_{0}={n}^{2}$ ($n=$ 0, 1, 2, …), where the effects of the fractionalorder derivative on the system damping and stiffness are formulated as equivalent linear damping coefficient (ELDC) and equivalent linear stiffness coefficient (ELSC). In section 3, the comparisons between the boundaries obtained by the approximate analytical solution and the numerical method verify the correctness and satisfactory precision of the analytical solution. The following analysis is focused on the effects of the fractional parameters on the stability boundaries in the $\delta \epsilon $ plane. At last the detail results are summarized and the conclusions are made.
2. Fractionalorder Mathieu equation
In this paper, we shall consider the fractionalorder Mathieu equation as:
where $2\zeta $, $\delta $, and $2\epsilon \mathrm{c}\mathrm{o}\mathrm{s}2t$ are the system linear damping coefficient, constant stiffness coefficient and the periodic timevarying coefficient respectively. ${D}^{p}\left[x\right(t\left)\right]$ is the $p$order derivative of $x\left(t\right)$ to $t$ with the fractional coefficient ${K}_{1}$ (${K}_{1}>0$) and the fractional order $p$ ($0\le p\le 1$). There are several definitions for fractionalorder derivative, such as GrünwaldLetnikov, RiemannLiouville and Caputo definitions [17]. Here Caputo’s definition is adopted with the form as:
where $\mathrm{\Gamma}\left(y\right)$ is Gamma function satisfying $\mathrm{\Gamma}\left(y+1\right)=y\mathrm{\Gamma}\left(y\right)$.
3. Approximately analytical solution and the stability boundaries
The approximately analytical solution and the stability boundaries are determined under the conditions of small damping coefficient. Using the parameters transformation:
Eq. (1) becomes:
Introducing the fast and slow time scales as follows:
and the differentiation operators:
the ordinary derivative can be transformed into the partial derivatives [34] as follows:
$\frac{{d}^{2}}{d{t}^{2}}={D}_{0}^{2}+2\epsilon {D}_{0}{D}_{1}+\dots .$
The approximate solution of Eq. (3) can be written as:
Applying the modified LP method, we introduce the following variable transformation [3435]:
Substituting Eq. (5), Eq. (6a), and Eq. (6b) into Eq. (3) and separating the terms in the yielded equation based on the power of $\epsilon $, we derive the following equations:
$2\mu {D}_{1}{x}_{0}\delta {}_{2}{}^{}{x}_{0}\delta {}_{1}{}^{}{x}_{1}2\mathrm{c}\mathrm{o}\mathrm{s}2T{}_{0}{}^{}{x}_{1}k{D}_{0}^{p}{x}_{1}k{D}_{1}^{p}{x}_{0}.$
The solution of Eq. (7) can be written as:
where $cc$ denotes the complex conjugate of the preceding term. The approximately analytical solution and the boundary will be presented in the following parts according to the different constant stiffness coefficient ${\delta}_{0}$.
3.1. ${\mathit{\delta}}_{0}=0$
According to Eq. (7) and Eq. (10), the solution can be written as:
where $a$ is an arbitrary constant and determined by the initial condition. Substituting Eq. (11) into Eq. (8) one could obtain:
In order to eliminate the secular terms in ${x}_{1}$, one must put:
The solution of Eq. (8) is:
Based on the formula for fractional derivative [36], ${D}^{p}{e}^{i\lambda t}$ could be approximately reduced to:
According to Eq. (13) and the Euler formula, one can obtain:
Substituting Eq. (11), Eq. (12), Eq. (13) and Eq. (14) into Eq. (9), we obtain:
To eliminate the secular terms, we arrive at:
and:
The solution of Eq. (9) is derived as:
where the formula:
is used in Eq. (16b). One can get another form of Eq. (16b):
Substituting Eq. (12a) and Eq. (16a) into Eq. (6b), the transition curve of the stability boundary near to ${\delta}_{0}=0$ is given by:
Substituting Eq. (11), Eq. (12b) and Eq. (18) into Eq. (6a), one could obtain the corresponding periodic solution with period $\pi $ on this curve:
$\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}+\frac{{\mathrm{\epsilon}}^{2}a}{32}\mathrm{c}\mathrm{o}\mathrm{s}4t+O\left({\mathrm{\epsilon}}^{3}\right).$
From Eq. (19), we can know that the result of the transition curve is independent of the fractional coefficient ${K}_{1}$ and the fractional order $p$. In order to analyze the effects of the fractionalorder derivative, we must obtain the thirdorder approximate solution:
$\delta {}_{3}{}^{}{x}_{0}\delta {}_{2}{}^{}{x}_{1}\delta {}_{1}{}^{}{x}_{2}2\mu {D}_{0}{x}_{2}2\mu {D}_{1}{x}_{1}2\mu {D}_{2}{x}_{0}$
$k{D}_{0}^{p}{x}_{2}k{D}_{1}^{p}{x}_{1}k{D}_{2}^{p}{x}_{0}.$
As the similar procedure, we arrive at:
and:
based on the eliminating procedure for the secular term.
In the case ${\delta}_{0}=0$, the thirdorder approximate solution could be established as:
where ${\mathrm{\Delta}}_{0}=\delta +{\epsilon}^{2}\frac{{2}^{p}{K}_{1}}{8}\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{p\pi}{2}\right)$, which is defined as the equivalent linear stiffness coefficient (ELSC) in the case ${\delta}_{0}=0$.
3.2. ${\mathit{\delta}}_{0}=\text{1}$
In the case ${\delta}_{0}=\text{1}$, the solution of Eq. (7) can be written as:
Using Eq. (13), one can get:
Substituting Eq. (24) and Eq. (25) into Eq. (8), we can obtain:
After eliminating the secular terms, we arrive at:
where ${D}_{1}A$ stands for the derivative with the respect to the slow time scale ${T}_{1}$, and $\stackrel{}{A}$ means the conjugate of $A$. In this paper, we want to obtain the periodic solutions, so that the hypothesis ${D}_{1}A=$0 is reasonable.
Substituting $A=a/2b/2i$ and Eq. (17) into Eq. (27), and after separating the real and imaginary parts of Eq. (27), we obtain:
The necessary conditions for existing nonzero solution about $(a,b)$ in Eq. (28) is:
That is:
Then, we can get:
In this case, the particular solution of Eq. (8) is:
Substituting Eq. (24) and Eq. (32) into Eq. (9), we can obtain:
$+\left[\frac{3}{4}i{D}_{1}A\frac{3}{4}i\mu A\frac{{\delta}_{1}A}{8}\frac{(3i{)}^{p}kA}{8}\right]{e}^{i3{T}_{0}}\frac{A}{8}{e}^{i5{T}_{0}}+cc.$
After eliminating the secular terms, one can get:
As the similar approach, we can suppose ${D}_{1}A=0$, ${D}_{2}A=0$ and ${D}_{1}^{p}A=0$. From Eq. (34) we can obtain:
In this case, the particular solution of Eq. (9) is:
Substituting Eq. (31) and Eq. (35) into Eq. (6b), the two transition curves emanating from ${\delta}_{0}=\text{1}$ are given by:
where ${\mathrm{\Delta}}_{1}=\delta +{K}_{1}\mathrm{c}\mathrm{o}\mathrm{s}\left(p\pi /2\right)$ and ${C}_{1}=2\zeta +{K}_{1}\mathrm{s}\mathrm{i}\mathrm{n}\left(p\pi /2\right)$. Here ${\mathrm{\Delta}}_{1}$ and ${C}_{1}$ are defined as the equivalent linear stiffness coefficient and the equivalent linear damping coefficient respectively in the case ${\delta}_{0}=\text{1}$.
Substituting Eq. (24), Eq. (32) and Eq. (36) into Eq. (6a), one could obtain the corresponding periodic solution with period 2$\pi $ on these curves:
Substituting Eq. (17), Eq. (37) and $A=a/2b/2i$ into Eq. (38), one can get another form of Eq. (38) with the original parameters:
$+\frac{{3}^{p}{K}_{1}\epsilon}{64}\left\{\left[\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{p\pi}{2}\right)a+\mathrm{s}\mathrm{i}\mathrm{n}\left(\frac{p\pi}{2}\right)b\right]\mathrm{c}\mathrm{o}\mathrm{s}3t+\left[\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{p\pi}{2}\right)b\mathrm{s}\mathrm{i}\mathrm{n}\left(\frac{p\pi}{2}\right)a\right]\mathrm{s}\mathrm{i}\mathrm{n}3t\right\}$
$+\frac{\epsilon}{64}\left[\pm \sqrt{{\epsilon}^{2}{{C}_{1}}^{2}}{\mathrm{\Delta}}_{1}\right]\left(a\mathrm{c}\mathrm{o}\mathrm{s}3t+b\mathrm{s}\mathrm{i}\mathrm{n}3t\right)+\frac{{\epsilon}^{2}}{192}\left(a\mathrm{c}\mathrm{o}\mathrm{s}5t+b\mathrm{s}\mathrm{i}\mathrm{n}5t\right)+O\left({\epsilon}^{3}\right).$
3.3. ${\mathit{\delta}}_{0}=\text{4}$
In order to analyze the effect of damping ratio near $\delta =\text{4}$, one should assumed $\zeta =O\left({\epsilon}^{2}\right)$. Using the transformation:
Eq. (1) becomes:
Eq. (8) and Eq. (9) become:
$\delta {}_{2}{}^{}{x}_{0}\delta {}_{1}{}^{}{x}_{1}2\mathrm{c}\mathrm{o}\mathrm{s}2{T}_{0}{x}_{1}.$
In the case ${\delta}_{0}=\text{4}$, the solution of Eq. (7) can be written as:
Substituting Eq. (43) into Eq. (41), we get:
To eliminate the secular terms, we arrive at $4i{D}_{1}A{\delta}_{1}A=0$. Due to ${D}_{1}^{p}A=0$, one can get ${\delta}_{1}=0$. The solution of Eq. (41) is:
Substituting Eq. (43) and Eq. (45) into Eq. (42), we get:
$\frac{2i{D}_{1}A}{3}{e}^{i4{T}_{0}}\frac{A}{12}{e}^{i6{T}_{0}}+cc.$
To eliminate the secular terms, we arrive at:
As the same procedure in the case ${\delta}_{0}=\text{1}$, we get:
That will result into:
Substituting the above results into Eq. (6b), the two transition curves emanating from ${\delta}_{0}=\text{4}$ are given by:
where:
${\mathrm{\Delta}}_{2}$, ${C}_{2}$ are defined as the equivalent linear stiffness coefficient (ELSC) and the equivalent linear damping coefficient (ELDC) respectively in the case ${\delta}_{0}=\text{4}$.
Substituting the above results into Eq. (6a), one could establish the corresponding periodic solution with period $\pi $ on these curves:
Transforming the exponent form into trigonometric one, we can obtain:
By analyzing the ELSC and the ELDC in three cases near ${\delta}_{0}={n}^{2}$ ($n=$ 0, 1, 2), the general forms of the ELSC and the ELDC are given by:
Here we can analyze some peculiar cases. If taking ${K}_{1}=0$, Eq. (1) becomes the classic Mathieu equation. The stability boundaries and the corresponding periodic solutions near ${\delta}_{0}={n}^{2}$ ($n=$ 0, 1, 2) are identical with the results in the classic monograph [3435]. If taking $p=\text{1}$ or $p=0$, the fractionalorder derivative will become the linear damping or the linear stiffness. The stability boundaries and the corresponding periodic solutions are identical with the results in the references [24] about the Mathieu equation. Therefore, it is indirectly verified that the correctness of this method and the results.
From the three cases, it can be concluded that the transition curve of the secondorder approximate solution near ${\delta}_{0}=\text{0}$ is independent of the fractionalorder derivative. But the fractional parameters will affect the transition curve of the thirdorder approximate solution near ${\delta}_{0}=\text{0}$ by the forms of the ELSC and the ELDC. It could also be found that the fractionalorder parameters will affect the transition curves of the secondorder approximate solution near ${\delta}_{0}=\text{1}$ and ${\delta}_{0}=\text{4}$ by the forms of the ELSC and the ELDC. In the three cases near ${\mathrm{\delta}}_{0}={n}^{2}$ ($n=$ 0, 1, 2), the ELSC and the ELDC could be arrived at the general forms shown in Eq. (52). Some parts of the corresponding periodic solutions on these boundaries can also be established by the forms of the ELSC and the ELDC.
From Eq. (52), it could be concluded that fractionalorder parameters ${K}_{1}$ and $p$ have important influence on the ELSC and the ELDC. It is easy to find that the equivalent linear damping and stiffness coefficients are all monotonically increasing function of the fractional coefficient ${K}_{1}$. The more important is the influence of the fractional order $p$ on the ELSC and the ELDC. When $p\to 0$, the ELDC will arrive at the minimum value as $2\zeta $, and the ELSC is $\delta +{K}_{1}$, so that the fractionalorder derivative is completely equivalent to the linear stiffness. On the contrary, when $p\to \text{1}$, the ELDC will be the maximum value as $2\zeta +{K}_{1}$, and the ELSC is $\delta $, so that the fractionalorder derivative is completely equivalent to the linear damping.
4. Numerical simulations
4.1. Comparison between the approximate analytical solution and the numerical results
In order to verify the validity of the method and the results in this paper, numerical results of the original equation are presented to compare the differences between the numerical simulations and the approximate analytical solutions as the next step.
In the $\delta \epsilon $ plane, we restrict the two parameters as $\epsilon \in $[0 0.4], and $\delta \in $[–0.2 4.5]. The sample step of $\epsilon $ and $\delta $ are selected as 0.01 and 0.005 respectively. Each point of the $\delta \epsilon $ plane is selected as the $\delta $ and $\epsilon $ value in Eq. (1) respectively, and used to calculate the numerical results. To determine the stability of each point, we should analyze the amplitude variation with a long time response. Then the stability boundaries near ${\delta}_{0}={n}^{2}$ ($n=$ 0, 1, 2, …) could be determined.
The numerical method is [17]:
where ${t}_{l}=lh$ is the time sample points, $h$ is the sample step, ${C}_{j}^{p}$ is the fractional binomial coefficient with the iterative relationship as:
Here we select $h=\mathrm{}$0.001 and the total computation time is 500 s.
An illustrative example system is studied herein as defined by the basic system parameters: $\zeta =$ 0.005, ${K}_{1}=$ 0.005, $p=$ 0.5. Based on the stability boundaries determined by Eq. (23), Eq. (37) and Eq. (48), one could analytically obtain those curves on the stability boundaries shown in Fig. 1, where the black solid lines are denoted for the approximate analytical solution. In Fig. 1, the red dots are for the unstable points and the blue dots are for the stable points, where the transition curves separating the red region from the blue region are the stability boundaries by numerical integration. From the observation of Fig. 1, we could conclude that the approximately analytical solution agrees very well with the numerical results and could present satisfactory precision.
Fig. 1The curves of the numerical results and the approximate analytical solution on the stability boundaries
4.2. Effects of fractionalorder parameters on the stability boundaries
In order to illustrate the effects of fractionalorder parameters on the stability boundaries, the fractional order $p$ are selected a set of values as 0, 0.2, 0.5, 0.8 and 1 respectively. The transition curves near ${\delta}_{0}=\text{1}$, ${\delta}_{0}=\text{4}$ and ${\delta}_{0}=\text{0}$ are shown in Fig. 2, Fig. 3 and Fig. 4 respectively. From the observation of the changes about the curves in Fig. 2 to Fig. 4, it is found that the increase of the fractional order $p$ could make the ELDC larger, which would result into the upwards moving of the stability boundaries. That is to say, the unstable region becomes smaller and the stable region becomes larger. Simultaneously, it could also be found that the increase of the fractional order $p$ would make the ELSC smaller, which could result into the rightwards moving of the stability boundaries. Therefore, it could be concluded that fractional order $p$ has important influence on the stability boundaries.
Fig. 2The effects of the fractional order p on the stability boundaries for δ0= 1 where ζ= 0.05 and K1= 0.05
Fig. 3The effects of the fractional order p on the stability boundaries for δ0= 4 where ζ= 0.005 and K1= 0.005
Fig. 4The effects of the fractional order p on the stability boundaries for δ0=0 where ζ= 0.5 and K1= 0.5
Fig. 5The effects of the fractional coefficient K1 on the stability boundaries for δ0= 1 where ζ= 0.05 and p= 0.5
Fig. 6The effects of the fractional coefficient K1 on the stability boundaries for δ0= 4 where ζ= 0.005 and p= 0.5
Fig. 7The effects of the fractional coefficient K1 on the stability boundaries for δ0=0 where p= 0.5
The fractional coefficient ${K}_{1}$ are selected some other different values, and the results are shown in Fig. 5 to Fig. 7. From the observation of those curves in Fig. 5 to Fig. 7, it is found that the increase of the fractional coefficient ${K}_{1}$ could make the ELSC larger, which would move the transition curves to the left. Simultaneously, it could be also found that the increase of the fractional coefficient ${K}_{1}$ will make the ELDC lager, which could move the transition curves to the upwards and make the unstable region smaller. Therefore, it could be concluded that the fractional coefficient ${K}_{1}$ has also important influence on the stability boundaries.
5. Conclusions
The dynamical characteristics of Mathieu equation with fractionalorder derivative is studied by the LindstedtPoincare method and the multiplescale method. The stability boundaries and the corresponding periodic solutions on these boundaries for the constant stiffness ${\delta}_{0}={n}^{2}$ ($n=$ 0, 1, 2, …), are analytically obtained. The effects of the fractionalorder parameters on the stability boundaries and the corresponding periodic solutions, including the fractional coefficient and the fractional order, are characterized by the equivalent linear damping coefficient (ELDC) and the equivalent linear stiffness coefficient (ELSC). The comparisons between the transition curves on the boundaries obtained by the approximate analytical solution and the numerical method verify the correctness and satisfactory precision of the analytical solution. It has been illustrated that the different fractionalorder parameters have important effects on the stability boundaries. These results are very helpful to design, analyze or control this kind of system, and could present beneficial reference to the similar fractionalorder system.
References

Miller K. S., Ross B. An Introduction to the Fractional Calculus and Fractional Differential Equations. Wiley, New York, 1993.

Podlubny I. Fractional Differential Equations. Academic, London, 1998.

Kilbas A. A., Srivastava H. M., Trujillo J. J. Theory and Applications of Fractional Differential Equations. Elsevier, Amsterdam, 2006.

Das S. Functional Fractional Calculus for System Identification and Controls. SpringerVerlag, Berlin, 2008.

Caponetto R., Dongola G., Fortuna L., Petras I. Fractional Order Systems: Modeling and Control Applications. World Scientific, New Jersey, 2010.

Monje C. A., Chen Y. Q., Vinagre B. M., Xue D. Y., Feliu V. Fractionalorder Systems and Controls: Fundamentals and Applications. SpringerVerlag, London, 2010.

Petras I. Fractionalorder Nonlinear Systems: Modeling, Analysis and Simulation. Higher Education Press, Beijing, 2011.

Shen Y. J., Yang S. P., Xing H. J., Gao G. S. Primary resonance of Duffing oscillator with fractionalorder derivative. Communications in Nonlinear Science and Numerical Simulation, Vol. 17, Issue 7, 2012, p. 30923100.

Shen Y. J., Yang S. P., Xing H. J., et al. Primary resonance of Duffing oscillator with two kinds of fractionalorder derivatives. International Journal of NonLinear Mechanics, Vol. 47, Issue 9, 2012, p. 975983.

Shen Y. J., Wei P., Yang S. P. Primary resonance of fractionalorder van der Pol oscillator. Nonlinear Dynamics, Vol. 77, Issue 4, 2014, p. 16291642.

Shen Y. J., Yang S. P., Sui C. Y. Analysis on limit cycle of fractionalorder van der Pol oscillator. Chaos, Solitons & Fractals, Vol. 67, 2014, p. 94102.

Gorenflo R., AbdelRehim E. A. Convergence of the GrünwaldLetnikov scheme for timefractional diffusion. Journal of Computational and Applied Mathematics, Vol. 205, Issue 2, 2007, p. 871881.

Jumarie G. Modified RiemannLiouville derivative and fractional Taylor series of non differentiable functions further results. Computers and Mathematics with Applications, Vol. 51, Issue 910, 2006, p. 13671376.

Ishteva M., Scherer R., Boyadjiev L. On the Caputo operator of fractional calculus and CLaguerre functions. Mathematical Sciences Research Journal, Vol. 9, Issue 6, 2005, p. 161170.

Agnieszka B. Malinowska, Delfim F. M. Torres Fractional calculus of variations for a combined Caputo derivative. Fractional Calculus and Applied Analysis, Vol. 14, Issue 4, 2011, p. 523537.

Chen L. C., Zhu W. Q. Stochastic dynamics and fractional optimal control of quasiintegrable Hamiltonian systems with fractional derivative damping. Fractional Calculus and Applied Analysis, Vol. 16, Issue 1, 2013, p. 189225.

Chen L. C., Zhu W. Q. The first passage failure of SDOF strongly nonlinear stochastic system with fractional derivative damping. Journal of Vibration Control, Vol. 15, Issue 8, 2009, p. 12471266.

Chen L. C., Zhu W. Q. Stochastic jump and bifurcation of Duffing oscillator with fractional derivative damping under combined harmonic and white noise excitations. International Journal of Nonlinear Mechanics, Vol. 46, Issue 12, 2011, p. 13241329.

Chen L. C., Li H. F., Li Z. S., Zhu W. Q. First passage failure of singledegreeoffreedom nonlinear oscillators with fractional derivative. Journal of Vibration and Control, Vol. 19, Issue 14, 2013, p. 21542163.

Chen L. C., Wang W. H., Li Z. S., Zhu W. Q. Stationary response of Duffing oscillator with hardening stiffness and fractional derivative. International Journal of Nonlinear Mechanics, Vol. 48, Issue 1, 2013, p. 4450.

Wang Z. H., Hu H. Y. Stability of a linear oscillator with damping force of fractionalorder derivative. Science in China G: Physics, Mechanics and Astronomy, Vol. 39, Issue 10, 2009, p. 14951502.

Wang Z. H., Du M. L. Asymptotical behavior of the solution of a SDOF linear fractionally damped vibration system. Shock and Vibration, Vol. 18, Issues 12, 2011, p. 257268.

Li C. P., Deng W. H. Remarks on fractional derivatives. Applied Mathematics and Computation, Vol. 187, Issue 2, 2007, p. 777784.

Cao J. X., Ding H. F., Li C. P. Implicit difference schemes for fractional diffusion equations. Communication on Applied Mathematics and Computation, Vol. 27, Issue 1, 2013, p. 6174.

Zeng F. H., Li C. P. Highorder finite difference methods for timefractional subdiffusion equation. Chinese Journal of Computational Physics, Vol. 30, Issue 4, 2013, p. 491500.

Chen A., Li C. P. Numerical algorithm for fractional calculus based on Chebyshev polynomial approximation. Journal of Shanghai Jiaotong University, Vol. 18, Issue 1, 2012, p. 4853.

Wahi P., Chatterjee A. Averaging oscillations with small fractional damping and delayed terms. Nonlinear Dynamics, Vol. 38, Issue 14, 2004, p. 322.

Xu Y., Li Y. G., Liu D., Jia W. T., Huang H. Responses of Duffing oscillator with fractional damping and random phase. Nonlinear Dynamics, Vol. 74, Issue 3, 2013, p. 745753.

Mclachlan N. W. Theory and Application of Mathieu Functions. Oxford University Press, London, 1951.

Ge Z. M., Yi C. X. Chaos in a nonlinear damped Mathieu system, in a nano resonator system and in its fractional order systems. Chaos, Solitons & Fractals, Vol. 32, Issue 1, 2007, p. 4261.

Ebaid A., ElSayed D. M. M., Aljoufi M. D. Fractional calculus model for damped Mathieu equation: approximate analytical solution. Applied Mathematical Sciences, Vol. 6, Issue 82, 2012, p. 40754080.

Rand R. H., Sah S. M., Suchorsky M. K. Fractional Mathieu equation. Communications in Nonlinear Science and Numerical Simulation, Vol. 15, Issue 11, 2010, p. 32543262.

Leung A. Y. T., Guo Z. J., Yang H. X. Transition curves and bifurcations of a class of fractional Mathieutype equations. International Journal of Bifurcation and Chaos, Vol. 22, 2012, p. 113.

Nayfeh A. H. Nonlinear Oscillations. Wiley, New York, 1979.

Nayfeh A. H. Introduction to Perturbation Techniques. Wiley, New York, 1981.

Rossikhin Y. A., Shitikova M. V. On fallacies in the decision between the Caputo and RiemannLiouville fractional derivative for the analysis of the dynamic response of a nonlinear viscoelastic oscillator. Mechanics Research Communications, Vol. 45, 2012, p. 2227.
About this article
The authors are grateful to the support by National Natural Science Foundation of China (No. 11372198), the Program for New Century Excellent Talents in University (NCET110936), the Cultivation Plan for Innovation Team and Leading Talent in Colleges and Universities of Hebei Province (LJRC018), the Program for Advanced Talent in the Universities of Hebei Province (GCC2014053), and the Program for Advanced Talent in Hebei (A201401001).