Dynamical analysis of fractional-order Mathieu equation

Wen, Shaofang; Shen, Yongjun; Li, Xianghong; Yang, Shaopu; Xing, Haijun

Journal of Vibroengineering

Browse Journal

Submit article

Published: 15 August 2015

Check for updates

Dynamical analysis of fractional-order Mathieu equation

Shaofang Wen¹

Yongjun Shen²

Xianghong Li³

Shaopu Yang⁴

Haijun Xing⁵

¹Transportation Institute, Shijiazhuang Tiedao University, Shijiazhuang, China

^{2, 4, 5}Department of Mechanical Engineering, Shijiazhuang Tiedao University, Shijiazhuang, China

³Department of Mathematics and Physics, Shijiazhuang Tiedao University, Shijiazhuang, China

Corresponding Author:

Yongjun Shen

Cite the article Download PDF

Downloads 1322

Abstract

The dynamical characteristics of Mathieu equation with fractional-order derivative is analytically studied by the Lindstedt-Poincare method and the multiple-scale method. The stability boundaries and the corresponding periodic solutions on these boundaries for the constant stiffness $δ_{0} = n^{2}$ ( $n$ = 0, 1, 2, …), are analytically obtained. The effects of the fractional-order 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 $δ - ε$ 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 fractional-order system.

1. Introduction

In 1695, the French mathematician Hospital and German mathematician Leibniz put forward the concept of fractional-order calculus for the first time. In the following centuries, the theory about fractional-order calculus including the definition, characteristics, calculation of fractional-order calculus and the relationship with integer-order calculus, etc, had been developing rapidly and a series of meaningful results had been obtained [1-7].

In recent years, fractional-order 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. [8-11] investigated several linear and nonlinear fractional-order oscillators by averaging method, and found that the fractional-order 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 fractional-order calculus for Grünwald-Letnikov, Riemann-Liouville and Caputo. Chen and Zhu [16] reviewed the analytical methods and control strategies for quasi-integrable Hamiltonian systems with fractional-order derivatives, and pointed out some possible developing issues. They [17-20] also studied some nonlinear fractional-order system with different kinds of noise, and obtained some important statistic properties of the fractional-order system. Wang et al. [21-22] investigated a linear single degree-of-freedom oscillator with fractional-order 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. [24-26] had done a lot of researches in the mathematical theory of fractional-order calculus, and also established some efficient numerical algorithms. Wahi and Chatterjee [27] studied a special linear single degree-of-freedom oscillator with fractional-order derivative by averaging method, and analyzed the effects of the fractional-orders derivative. Xu et al. [28] used Lindstedt-Poincaré method and the multiple-scale approach to investigate fractional-order 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 integer-order system, the fractional-order system is much closer to the real nature of the world, and has more advantages, such as strong ability of anti-noise, good robustness, high control precision and so on. Accordingly fractional-order calculus has the significant value for the engineering field.

The well-known 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 fractional-order system could generate different dynamical properties from the integer-order 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 fractional-order 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 fractional-order Mathieu equation had been studied by the method of harmonic balance [32-33], and Leung et al. also investigated a general version of fractional-order Mathieu-Duffing equation by harmonic balance method.

In this paper, the fractional-order Mathieu equation with linear damping is studied by the Lindstedt-Poincaré (L-P) method coupled with the multiple-scale method, and the effects of the fractional parameters on the dynamical properties of the fractional-order 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 $δ_{0} = n^{2}$ ( $n =$ 0, 1, 2, …), where the effects of the fractional-order 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 $δ - ε$ plane. At last the detail results are summarized and the conclusions are made.

2. Fractional-order Mathieu equation

In this paper, we shall consider the fractional-order Mathieu equation as:

1

\ddot{x} (t) + 2 ζ \dot{x} (t) + (δ + 2 ε c o s 2 t) x (t) + K_{1} D^{p} [x (t)] = 0, |ε| ≪ 1, 0 < ζ ≪ 1,

where $2 ζ$ , $δ$ , and $2 ε c o s 2 t$ are the system linear damping coefficient, constant stiffness coefficient and the periodic time-varying coefficient respectively. $D^{p} [x (t)]$ is the $p$ -order derivative of $x (t)$ to $t$ with the fractional coefficient $K_{1}$ ( $K_{1} > 0$ ) and the fractional order $p$ ( $0 \leq p \leq 1$ ). There are several definitions for fractional-order derivative, such as Grünwald-Letnikov, Riemann-Liouville and Caputo definitions [1-7]. Here Caputo’s definition is adopted with the form as:

2

D^{p} [x (t)] = \frac{1}{Γ (1 - p)} \int_{0}^{t} \frac{x^{'} (u)}{{(t - u)}^{p}} d u,

where $Γ (y)$ is Gamma function satisfying $Γ (y + 1) = y Γ (y)$ .

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:

ζ = ε μ, μ = O (1), K_{1} = ε k, k = O (1) .

Eq. (1) becomes:

3

\ddot{x} (t) + 2 ε μ \dot{x} (t) + (δ + 2 ε c o s 2 t) x (t) + ε k D^{p} [x (t)] = 0 .

Introducing the fast and slow time scales as follows:

T_{r} \overset{d e f}{=} ε^{r} t, r = 0, 1, 2, \dots,

and the differentiation operators:

4

D_{0} = \frac{\partial}{\partial T_{0}}, D_{1} = \frac{\partial}{\partial T_{1}}, D_{2} = \frac{\partial}{\partial T_{2}}

the ordinary derivative can be transformed into the partial derivatives [34] as follows:

5

\begin{matrix} \frac{d}{d t} = D_{0} + ε D_{1} + \dots, \end{matrix}

\frac{d^{2}}{d t^{2}} = D_{0}^{2} + 2 ε D_{0} D_{1} + \dots .

The approximate solution of Eq. (3) can be written as:

6a

x = x_{0} (T_{0}, T_{1}, T_{2}) + ε x_{1} (T_{0}, T_{1}, T_{2}) + ε^{2} x_{2} (T_{0}, T_{1}, T_{2}) + \dots .

Applying the modified L-P method, we introduce the following variable transformation [34-35]:

6b

δ = δ_{0} + ε δ_{1} + ε^{2} δ_{2} + ε^{3} δ_{3} + \dots .

Substituting Eq. (5), Eq. (6a), and Eq. (6b) into Eq. (3) and separating the terms in the yielded equation based on the power of $ε$ , we derive the following equations:

7

O (ε^{0}) : D_{0}^{2} x_{0} + δ {}_{0} x_{0} = 0,

8

O (ε^{1}) : D_{0}^{2} x_{1} + δ {}_{0} x_{1} = - 2 D_{0} D_{1} x_{0} - 2 μ D_{0} x_{0} - δ {}_{1} x_{0} - 2 c o s 2 T {}_{0} x_{0} - k D_{0}^{p} x_{0},

9

O (ε^{2}) : D_{0}^{2} x_{2} + δ {}_{0} x_{2} = - 2 D_{0} D_{2} x_{0} - 2 D_{0} D_{1} x_{1} - D_{1}^{2} x_{0} - 2 μ D_{0} x_{1}

- 2 μ D_{1} x_{0} - δ {}_{2} x_{0} - δ {}_{1} x_{1} - 2 c o s 2 T {}_{0} x_{1} - k D_{0}^{p} x_{1} - k D_{1}^{p} x_{0} .

The solution of Eq. (7) can be written as:

10

x_{0} = A (T_{1}, T_{2}) e^{i \sqrt{δ_{0}} T_{0}} + c c,

where $c c$ 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 $δ_{0}$ .

3.1. $δ_{0} = 0$

According to Eq. (7) and Eq. (10), the solution can be written as:

11

x_{0} = a = c o n s t,

where $a$ is an arbitrary constant and determined by the initial condition. Substituting Eq. (11) into Eq. (8) one could obtain:

D_{0}^{2} x_{1} = - a (δ_{1} + 2 c o s 2 T_{0}) .

In order to eliminate the secular terms in $x_{1}$ , one must put:

12a

δ_{1} = 0 .

The solution of Eq. (8) is:

12b

x_{1} = \frac{a}{2} c o s 2 T_{0} .

Based on the formula for fractional derivative [36], $D^{p} e^{i λ t}$ could be approximately reduced to:

13

D^{p} e^{i λ t} = (i λ)^{p} e^{i λ t} .

According to Eq. (13) and the Euler formula, one can obtain:

14

k D_{0}^{p} x_{1} = k D_{0}^{p} (\frac{a}{2} c o s 2 T_{0}) = \frac{a}{2} k D_{0}^{p} (\frac{e^{i 2 T_{0}} + e^{- i 2 T_{0}}}{2}) = \frac{k a}{4} [(2 i)^{p} e^{i 2 T_{0}} + (- 2 i)^{p} e^{- i 2 T_{0}}] .

Substituting Eq. (11), Eq. (12), Eq. (13) and Eq. (14) into Eq. (9), we obtain:

15

D_{0}^{2} x_{2} = - a (δ_{2} + \frac{1}{2}) + 2 μ a s i n 2 T_{0} - \frac{a}{2} c o s 4 T_{0} - \frac{1}{4} k a [(2 i)^{p} e^{i 2 T_{0}} + (- 2 i)^{p} e^{- i 2 T_{0}}] .

To eliminate the secular terms, we arrive at:

16

- a (δ_{2} + \frac{1}{2}) = 0,

and:

16a

δ_{2} = - \frac{1}{2} .

The solution of Eq. (9) is derived as:

16b

x_{2} = - \frac{μ a}{2} s i n 2 T_{0} + \frac{a}{32} c o s 4 T_{0} + \frac{1}{16} k a (2 i)^{p} e^{i 2 T_{0}} + \frac{1}{16} k a (- 2 i)^{p} e^{- i 2 T_{0}},

where the formula:

17

i^{p} = (e^{\frac{i π}{2}})^{p} = e^{\frac{i p π}{2}} = c o s \frac{p π}{2} + i s i n \frac{p π}{2},

is used in Eq. (16b). One can get another form of Eq. (16b):

18

x {}_{2} = \frac{a}{2} [- μ s i n 2 T_{0} + \frac{2^{p} k}{4} c o s (\frac{p π}{2} + 2 T_{0})] + \frac{a}{32} c o s 4 T_{0} + O (ε^{3}) .

Substituting Eq. (12a) and Eq. (16a) into Eq. (6b), the transition curve of the stability boundary near to $δ_{0} = 0$ is given by:

19

δ (ε) = - \frac{ε^{2}}{2} + O (ε^{3}) .

Substituting Eq. (11), Eq. (12b) and Eq. (18) into Eq. (6a), one could obtain the corresponding periodic solution with period $π$ on this curve:

20

x = a + \frac{ε a}{2} \{[1 + \frac{2^{p} K {}_{1}}{4} c o s (\frac{p π}{2})] c o s 2 t - [ζ + \frac{2^{p} K {}_{1}}{4} s i n (\frac{p π}{2})] s i n 2 t\}

+ \frac{ε^{2} a}{32} c o s 4 t + O (ε^{3}) .

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 fractional-order derivative, we must obtain the third-order approximate solution:

21

O (ε^{3}) : D_{0}^{2} x_{3} + δ {}_{0} x_{3} = - 2 D_{1} D_{0} x_{2} - 2 D_{1} D_{2} x_{0} - D_{1}^{2} x_{1} - 2 D_{2} D_{0} x_{1} - 2 c o s 2 T {}_{0} x_{2}

- δ {}_{3} x_{0} - δ {}_{2} x_{1} - δ {}_{1} x_{2} - 2 μ D_{0} x_{2} - 2 μ D_{1} x_{1} - 2 μ 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:

- δ_{3} a - \frac{2^{p} k a}{8} c o s (\frac{p π}{2}) = 0

and:

22

δ_{3} = - \frac{2^{p} k}{8} c o s (\frac{p π}{2}),

based on the eliminating procedure for the secular term.

In the case $δ_{0} = 0$ , the third-order approximate solution could be established as:

23

Δ_{0} = - \frac{ε^{2}}{2} + O (ε^{4}),

where $Δ_{0} = δ + ε^{2} \frac{2^{p} K_{1}}{8} c o s (\frac{p π}{2})$ , which is defined as the equivalent linear stiffness coefficient (ELSC) in the case $δ_{0} = 0$ .

3.2. $δ_{0} = 1$

In the case $δ_{0} = 1$ , the solution of Eq. (7) can be written as:

24

x_{0} = A e^{i T_{0}} + c c .

Using Eq. (13), one can get:

25

k D_{0}^{p} x_{0} = k D_{0}^{p} (A e^{i T_{0}} + c c) = k i^{p} A e^{i T_{0}} + c c .

Substituting Eq. (24) and Eq. (25) into Eq. (8), we can obtain:

26

D_{0}^{2} x_{1} + x_{1} = (- 2 i D_{1} A - 2 i μ A - δ_{1} A - \bar{A} - k i^{p} A) e^{i T_{0}} - A e^{i 3 T_{0}} + c c .

After eliminating the secular terms, we arrive at:

27

- 2 i D_{1} A - 2 i μ A - δ_{1} A - \bar{A} - k i^{p} A = 0

where $D_{1} A$ stands for the derivative with the respect to the slow time scale $T_{1}$ , and $\bar{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 / 2 - b / 2 i$ and Eq. (17) into Eq. (27), and after separating the real and imaginary parts of Eq. (27), we obtain:

28

\{\begin{array}{l} [1 + δ_{1} + k c o s (\frac{p π}{2})] a + [2 μ + k s i n (\frac{p π}{2})] b = 0, \\ [2 μ + k s i n (\frac{p π}{2})] a + [1 - δ_{1} - k c o s (\frac{p π}{2})] b = 0 . \end{array}

The necessary conditions for existing nonzero solution about $(a, b)$ in Eq. (28) is:

29

d e t [\begin{matrix} 1 + δ_{1} + k c o s (\frac{p π}{2}) & 2 μ + k s i n (\frac{p π}{2}) \\ 2 μ + k s i n (\frac{p π}{2}) & 1 - δ_{1} - k c o s (\frac{p π}{2}) \end{matrix}] = 0 .

That is:

30

[1 + δ_{1} + k c o s (\frac{p π}{2})] [1 - δ_{1} - k c o s (\frac{p π}{2})] - {[2 μ + k s i n (\frac{p π}{2})]}^{2} = 0 .

Then, we can get:

31

δ_{1} = \pm \sqrt{1 - {[2 μ + k s i n (\frac{p π}{2})]}^{2}} - k c o s (\frac{p π}{2}) .

In this case, the particular solution of Eq. (8) is:

32

x_{1} = \frac{A}{8} e^{i 3 T_{0}} + c c .

Substituting Eq. (24) and Eq. (32) into Eq. (9), we can obtain:

33

D_{0}^{2} x_{2} + x_{2} = (- 2 i D_{2} A - D_{1}^{2} A - 2 μ D_{1} A - δ_{2} A - \frac{A}{8} - k D_{1}^{p} A) e^{i T_{0}}

+ [- \frac{3}{4} i D_{1} A - \frac{3}{4} i μ A - \frac{δ_{1} A}{8} - \frac{(3 i)^{p} k A}{8}] e^{i 3 T_{0}} - \frac{A}{8} e^{i 5 T_{0}} + c c .

After eliminating the secular terms, one can get:

34

- 2 i D_{2} A - D_{1}^{2} A - 2 μ D_{1} A - δ_{2} A - \frac{A}{8} - k D_{1}^{p} A = 0 .

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:

35

δ_{2} = - \frac{1}{8} .

In this case, the particular solution of Eq. (9) is:

36

x_{2} = (i \frac{3 μ A}{32} + \frac{δ_{1} A}{64} + i^{p} \frac{3^{p} k A}{64}) e^{i 3 T_{0}} + \frac{A}{192} e^{i 5 T_{0}} + c c .

Substituting Eq. (31) and Eq. (35) into Eq. (6b), the two transition curves emanating from $δ_{0} = 1$ are given by:

37

Δ_{1} = 1 \pm \sqrt{ε^{2} - C_{1}^{2}} - \frac{1}{8} ε^{2} + O (ε^{3}),

where $Δ_{1} = δ + K_{1} c o s (p π / 2)$ and $C_{1} = 2 ζ + K_{1} s i n (p π / 2)$ . Here $Δ_{1}$ and $C_{1}$ are defined as the equivalent linear stiffness coefficient and the equivalent linear damping coefficient respectively in the case $δ_{0} = 1$ .

Substituting Eq. (24), Eq. (32) and Eq. (36) into Eq. (6a), one could obtain the corresponding periodic solution with period 2 $π$ on these curves:

38

x = A e^{i T_{0}} + ε \frac{A}{8} e^{i 3 T_{0}} + ε^{2} (i \frac{3 μ A}{32} + \frac{δ_{1} A}{64} + i^{p} \frac{3^{p} k A}{64}) e^{i 3 T_{0}} + ε^{2} \frac{A}{192} e^{i 5 T_{0}} + c c + O (ε^{3}) .

Substituting Eq. (17), Eq. (37) and $A = a / 2 - b / 2 i$ into Eq. (38), one can get another form of Eq. (38) with the original parameters:

39

x = a c o s t + b s i n t + \frac{ε}{8} \{[(1 + \frac{δ}{8}) a + \frac{3 ζ b}{4}] c o s 3 t + [(1 + \frac{δ}{8}) b - \frac{3 ζ a}{4}] s i n 3 t\}

+ \frac{3^{p} K_{1} ε}{64} \{[c o s (\frac{p π}{2}) a + s i n (\frac{p π}{2}) b] c o s 3 t + [c o s (\frac{p π}{2}) b - s i n (\frac{p π}{2}) a] s i n 3 t\}

+ \frac{ε}{64} [\pm \sqrt{ε^{2} - {C_{1}}^{2}} - Δ_{1}] (a c o s 3 t + b s i n 3 t) + \frac{ε^{2}}{192} (a c o s 5 t + b s i n 5 t) + O (ε^{3}) .

3.3. $δ_{0} = 4$

In order to analyze the effect of damping ratio near $δ = 4$ , one should assumed $ζ = O (ε^{2})$ . Using the transformation:

ζ = ε^{2} \hat{μ}, \hat{μ} = O (1), K_{1} = ε^{2} \hat{k}, \hat{k} = O (1) .

Eq. (1) becomes:

40

\ddot{x} (t) + 2 ε^{2} \hat{μ} \dot{x} (t) + (δ + 2 ε c o s 2 t) x (t) + ε^{2} \hat{k} D^{p} [x (t)] = 0 .

Eq. (8) and Eq. (9) become:

41

O (ε^{1}) : D_{0}^{2} x_{1} + δ_{0} x_{1} = - 2 D_{0} D_{1} x_{0} - δ_{1} x_{0} - 2 c o s 2 T_{0} x_{0},

42

O (ε^{2}) : D_{0}^{2} x_{2} + δ {}_{0} x_{2} = - 2 D_{0} D_{2} x_{0} - 2 D_{0} D_{1} x_{1} - D_{1}^{2} x_{0} - 2 \hat{μ} D_{0} x_{0} - \hat{k} D_{0}^{p} x_{0}

- δ {}_{2} x_{0} - δ {}_{1} x_{1} - 2 c o s 2 T_{0} x_{1} .

In the case $δ_{0} = 4$ , the solution of Eq. (7) can be written as:

43

x_{0} = A e^{i 2 T_{0}} + c c .

Substituting Eq. (43) into Eq. (41), we get:

44

D_{0}^{2} x_{1} + 4 x_{1} = (- 4 i D_{1} A - δ_{1} A) e^{i 2 T_{0}} - A e^{i 4 T_{0}} + c c .

To eliminate the secular terms, we arrive at $- 4 i D_{1} A - δ_{1} A = 0$ . Due to $D_{1}^{p} A = 0$ , one can get $δ_{1} = 0$ . The solution of Eq. (41) is:

45

x_{1} = \frac{A}{12} e^{i 4 T_{0}} - \frac{A}{4} + c c .

Substituting Eq. (43) and Eq. (45) into Eq. (42), we get:

46

D_{0}^{2} x_{2} + 4 x_{2} = (- 4 i D_{2} A - D_{1}^{2} A - 4 i \hat{μ} A - 2^{p} i^{p} \hat{k} A - δ_{2} A + \frac{A}{6} - \frac{\bar{A}}{4}) e^{i 2 T_{0}}

- \frac{2 i D_{1} A}{3} e^{i 4 T_{0}} - \frac{A}{12} e^{i 6 T_{0}} + c c .

To eliminate the secular terms, we arrive at:

- 4 i D_{2} A - D_{1}^{2} A - 4 i \hat{μ} A - 2^{p} i^{p} \hat{k} A - δ_{2} A + \frac{A}{6} - \frac{\bar{A}}{4} = 0 .

As the same procedure in the case $δ_{0} = 1$ , we get:

d e t [\begin{matrix} - [2^{p} \hat{k} c o s (\frac{p π}{2})] + δ_{2} + \frac{1}{12} & - [2^{p} \hat{k} s i n (\frac{p π}{2}) + 4 \hat{μ}] \\ - [2^{p} \hat{k} s i n (\frac{p π}{2}) + 4 \hat{μ}] & [2^{p} \hat{k} c o s (\frac{p π}{2}) + δ_{2} - \frac{5}{12}] \end{matrix}] = 0 .

That will result into:

47

δ_{2} = \frac{1}{6} - 2^{p} \hat{k} c o s (\frac{p π}{2}) \pm \sqrt{\frac{1}{16} - [2^{p} \hat{k} s i n (\frac{p π}{2}) + 4 \hat{μ}]^{2}} .

Substituting the above results into Eq. (6b), the two transition curves emanating from $δ_{0} = 4$ are given by:

48

Δ_{2} = 4 + \frac{ε^{2}}{6} \pm \sqrt{\frac{ε^{4}}{16} - 4 C_{2}^{2}} + O (ε^{3}),

where:

49

Δ_{2} = δ + 2^{p} K_{1} c o s (\frac{p π}{2}), C_{2} = 2 ζ + 2^{p - 1} K_{1} s i n (\frac{p π}{2}),

$Δ_{2}$ , $C_{2}$ are defined as the equivalent linear stiffness coefficient (ELSC) and the equivalent linear damping coefficient (ELDC) respectively in the case $δ_{0} = 4$ .

Substituting the above results into Eq. (6a), one could establish the corresponding periodic solution with period $π$ on these curves:

50

x = A e^{i 2 T_{0}} + ε (\frac{A}{12} e^{i 4 T_{0}} - \frac{A}{4}) + ε^{2} \frac{A}{384} e^{i 6 T_{0}} + c c + O (ε^{3}) .

Transforming the exponent form into trigonometric one, we can obtain:

51

x = a c o s 2 t + b s i n 2 t - ε (\frac{a}{4} - \frac{a}{12} c o s 4 t - \frac{b}{12} s i n 4 t) + \frac{ε^{2}}{384} (a c o s 6 t + b s i n 6 t) + O (ε^{3}) .

By analyzing the ELSC and the ELDC in three cases near $δ_{0} = n^{2}$ ( $n =$ 0, 1, 2), the general forms of the ELSC and the ELDC are given by:

52a

Δ = δ + {(\sqrt{δ_{0}})}^{p} K_{1} c o s (\frac{p π}{2}),

52b

C = 2 ζ + {(\sqrt{δ_{0}})}^{p - 1} K_{1} s i n (\frac{p π}{2}) .

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 $δ_{0} = n^{2}$ ( $n =$ 0, 1, 2) are identical with the results in the classic monograph [34-35]. If taking $p = 1$ or $p = 0$ , the fractional-order 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 [2-4] 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 second-order approximate solution near $δ_{0} = 0$ is independent of the fractional-order derivative. But the fractional parameters will affect the transition curve of the third-order approximate solution near $δ_{0} = 0$ by the forms of the ELSC and the ELDC. It could also be found that the fractional-order parameters will affect the transition curves of the second-order approximate solution near $δ_{0} = 1$ and $δ_{0} = 4$ by the forms of the ELSC and the ELDC. In the three cases near $δ_{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 fractional-order 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 ζ$ , and the ELSC is $δ + K_{1}$ , so that the fractional-order derivative is completely equivalent to the linear stiffness. On the contrary, when $p \to 1$ , the ELDC will be the maximum value as $2 ζ + K_{1}$ , and the ELSC is $δ$ , so that the fractional-order 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 $δ - ε$ plane, we restrict the two parameters as $ε \in$ [0 0.4], and $δ \in$ [–0.2 4.5]. The sample step of $ε$ and $δ$ are selected as 0.01 and 0.005 respectively. Each point of the $δ - ε$ plane is selected as the $δ$ and $ε$ 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 $δ_{0} = n^{2}$ ( $n =$ 0, 1, 2, …) could be determined.

The numerical method is [1-7]:

53

D^{p} [x (t_{l})] \approx h^{- p} \sum_{j = 0}^{l} C_{j}^{p} x (t_{l - j}),

where $t_{l} = l h$ is the time sample points, $h$ is the sample step, $C_{j}^{p}$ is the fractional binomial coefficient with the iterative relationship as:

C_{0}^{p} = 1, C_{j}^{p} = (1 - \frac{1 + p}{j}) C_{j - 1}^{p} .

Here we select $h =$ 0.001 and the total computation time is 500 s.

An illustrative example system is studied herein as defined by the basic system parameters: $ζ =$ 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 fractional-order parameters on the stability boundaries

In order to illustrate the effects of fractional-order 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 $δ_{0} = 1$ , $δ_{0} = 4$ and $δ_{0} = 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

$The 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

$The 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

$The 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

$The 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

$The 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 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 fractional-order derivative is studied by the Lindstedt-Poincare method and the multiple-scale method. The stability boundaries and the corresponding periodic solutions on these boundaries for the constant stiffness $δ_{0} = n^{2}$ ( $n =$ 0, 1, 2, …), are analytically obtained. The effects of the fractional-order 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 fractional-order 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 fractional-order system.

References

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

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

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

Search CrossRef
Das S. Functional Fractional Calculus for System Identification and Controls. Springer-Verlag, Berlin, 2008.

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

Search CrossRef
Monje C. A., Chen Y. Q., Vinagre B. M., Xue D. Y., Feliu V. Fractional-order Systems and Controls: Fundamentals and Applications. Springer-Verlag, London, 2010.

Search CrossRef
Petras I. Fractional-order Nonlinear Systems: Modeling, Analysis and Simulation. Higher Education Press, Beijing, 2011.

Search CrossRef
Shen Y. J., Yang S. P., Xing H. J., Gao G. S. Primary resonance of Duffing oscillator with fractional-order derivative. Communications in Nonlinear Science and Numerical Simulation, Vol. 17, Issue 7, 2012, p. 3092-3100.

Search CrossRef
Shen Y. J., Yang S. P., Xing H. J., et al. Primary resonance of Duffing oscillator with two kinds of fractional-order derivatives. International Journal of Non-Linear Mechanics, Vol. 47, Issue 9, 2012, p. 975-983.

Search CrossRef
Shen Y. J., Wei P., Yang S. P. Primary resonance of fractional-order van der Pol oscillator. Nonlinear Dynamics, Vol. 77, Issue 4, 2014, p. 1629-1642.

Search CrossRef
Shen Y. J., Yang S. P., Sui C. Y. Analysis on limit cycle of fractional-order van der Pol oscillator. Chaos, Solitons & Fractals, Vol. 67, 2014, p. 94-102.

Search CrossRef
Gorenflo R., Abdel-Rehim E. A. Convergence of the Grünwald-Letnikov scheme for time-fractional diffusion. Journal of Computational and Applied Mathematics, Vol. 205, Issue 2, 2007, p. 871-881.

Search CrossRef
Jumarie G. Modified Riemann-Liouville derivative and fractional Taylor series of non differentiable functions further results. Computers and Mathematics with Applications, Vol. 51, Issue 9-10, 2006, p. 1367-1376.

Search CrossRef
Ishteva M., Scherer R., Boyadjiev L. On the Caputo operator of fractional calculus and C-Laguerre functions. Mathematical Sciences Research Journal, Vol. 9, Issue 6, 2005, p. 161-170.

Search CrossRef
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. 523-537.

Search CrossRef
Chen L. C., Zhu W. Q. Stochastic dynamics and fractional optimal control of quasi-integrable Hamiltonian systems with fractional derivative damping. Fractional Calculus and Applied Analysis, Vol. 16, Issue 1, 2013, p. 189-225.

Search CrossRef
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. 1247-1266.

Search CrossRef
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 Non-linear Mechanics, Vol. 46, Issue 12, 2011, p. 1324-1329.

Search CrossRef
Chen L. C., Li H. F., Li Z. S., Zhu W. Q. First passage failure of single-degree-of-freedom nonlinear oscillators with fractional derivative. Journal of Vibration and Control, Vol. 19, Issue 14, 2013, p. 2154-2163.

Search CrossRef
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 Non-linear Mechanics, Vol. 48, Issue 1, 2013, p. 44-50.

Search CrossRef
Wang Z. H., Hu H. Y. Stability of a linear oscillator with damping force of fractional-order derivative. Science in China G: Physics, Mechanics and Astronomy, Vol. 39, Issue 10, 2009, p. 1495-1502.

Search CrossRef
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 1-2, 2011, p. 257-268.

Search CrossRef
Li C. P., Deng W. H. Remarks on fractional derivatives. Applied Mathematics and Computation, Vol. 187, Issue 2, 2007, p. 777-784.

Search CrossRef
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. 61-74.

Search CrossRef
Zeng F. H., Li C. P. High-order finite difference methods for time-fractional subdiffusion equation. Chinese Journal of Computational Physics, Vol. 30, Issue 4, 2013, p. 491-500.

Search CrossRef
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. 48-53.

Search CrossRef
Wahi P., Chatterjee A. Averaging oscillations with small fractional damping and delayed terms. Nonlinear Dynamics, Vol. 38, Issue 1-4, 2004, p. 3-22.

Search CrossRef
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. 745-753.

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

Search CrossRef
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. 42-61.

Search CrossRef
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. 4075-4080.

Search CrossRef
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. 3254-3262.

Search CrossRef
Leung A. Y. T., Guo Z. J., Yang H. X. Transition curves and bifurcations of a class of fractional Mathieu-type equations. International Journal of Bifurcation and Chaos, Vol. 22, 2012, p. 1-13.

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

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

Search CrossRef
Rossikhin Y. A., Shitikova M. V. On fallacies in the decision between the Caputo and Riemann-Liouville fractional derivative for the analysis of the dynamic response of a nonlinear viscoelastic oscillator. Mechanics Research Communications, Vol. 45, 2012, p. 22-27.

Search CrossRef

About this article

Received

01 February 2015

Accepted

11 April 2015

Published

15 August 2015

SUBJECTS

Chaos, nonlinear dynamics and applications

Keywords

fractional-order derivative

Mathieu equation

Lindstedt-Poincare method

multiple-scale method

stability boundaries

Acknowledgements

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 (NCET-11-0936), 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).

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.