# Hopf Bifurcations of a Stochastic Fractional-Order Van der Pol System

Authors: Xiaojun Liu; Ling Hong; Lixin Yang
Journal: Abstract and Applied Analysis
 (2014)

---

## Abstract

The Hopf bifurcation of a fractional-order Van der Pol (VDP for short) system with a random parameter is investigated. Firstly, the Chebyshev polynomial approximation is applied to study the stochastic fractional-order system. Based on the method, the stochastic system is reduced to the equivalent deterministic one, and then the responses of the stochastic system can be obtained by numerical methods. Then, according to the existence conditions of Hopf bifurcation, the critical parameter value of the bifurcation is obtained by theoretical analysis. Then, numerical simulations are carried out to verify the theoretical results.


---

## Body

## 1. Introduction

Fractional calculus is a topic of more than 300 years old. The idea of fractional calculus has been known since the regular calculus. However, its development is very rapid just in recent several decades due to its application in physics, engineering, secure communications, and so on. Meanwhile, it has been found that many fractional-order nonlinear systems can demonstrate chaotic behavior, such as fractional-order Chua circuit [1], fractional-order Lorenz system [2], and fractional-order Chen system [3]. These examples and many other similar samples perfectly clarify the importance of consideration and analysis of dynamical systems with fractional-order models [4–6].The VDP oscillator represents a nonlinear system with an interesting behavior that exhibits naturally in several applications. It has been used for study and design of many models including biological phenomena, such as the heartbeat, neurons, acoustic models, and radiation of mobile phones, and as a model of electrical oscillators. The deterministic VDP has very rich dynamical behaviors, such as saddle-node bifurcation, symmetry-breaking bifurcation, period-doubling bifurcation, Hopf bifurcation, and chaos. Recently, many researchers studied the bifurcation and periodic solutions of the VDP system in detail [7–11]. As research on fractional differential equations goes up, the VDP model with fractional order was proposed and investigated in [12, 13]. It was found that when the order is less than 1, the fractional-order VDP system also has rich dynamics similar to the corresponding integer-order one.In the real world, uncertainty due to measuring, materials, manufacturing, and assembling is inevitable. These uncertainties usually can be modeled as random parameters with certain statistics. A stochastic system means a system with random parameters of given statistics, which is a mathematical model close to the real world. There are several methods to deal with the dynamical problems of stochastic systems. The first one is the Monte Carlo method [14], which is simple and universal but usually involves a great deal of computational effort. The second one is the stochastic finite element method [15], which involves the least computation but is usually restricted to systems with random variables of small perturbations only. And the third one is orthogonal polynomials [16–20], which is based on the expansion theory of orthogonal polynomials, without the limitation of small perturbations. There are a lot of references about the analysis of stochastic systems via orthogonal polynomials. For example, in [21, 22], the Chebyshev polynomial approximation was applied to analyze the bifurcation and chaos of stochastic Duffing system and stochastic Duffing-Van der Pol system. The period-doubling bifurcation for double-well stochastic Duffing system was investigated via Lagurre polynomial approximation in [23]. Recently, the synchronization for a stochastic fractional-order system with a random parameter via Laguerre polynomial was analyzed in [24]. However, the literature on the investigation of a stochastic fractional-order VDP system is few.Inspired by the above discussion, in this paper, Hopf bifurcations of a fractional-order VDP system with a random parameter are studied. Firstly, the Chebyshev polynomial approximation method is applied to investigate the stochastic fractional-order system. Based on this method, the stochastic system is reduced to the equivalent deterministic one. Then, according to the existence conditions of Hopf bifurcation, the critical parameter value of the bifurcation is obtained by theoretical analysis. Meanwhile, numerical simulations are performed to verify the theoretical results. Besides, the Hopf bifurcation with the variation of derivative order is also observed by numerical computations.The paper is organized as follows. In Section2, the definitions for the fractional calculus and numerical algorithms are given. In Section 3, the Chebyshev polynomials are introduced as the basis of orthogonal polynomial approximation. In Section 4, the transformation of stochastic fractional-order VDP system into its equivalent deterministic one by Chebyshev polynomial approximation is shown. Section 5 is devoted to studying the Hopf bifurcation of the stochastic system. Finally, we summarize the results in Section 6.


## 2. Fundamentals of Fractional Derivative

### 2.1. Definition

Fractional calculus is a generalization of integration and differentiation to a noninteger-order integrodifferential operatorDtqa which is defined by
(1)Dtqa={dqdtqR(q)>01R(q)=0∫at(dτ)-qR(q)<0,
where q is the derivative order which can be a complex number and R(q) is the real part of q. The numbers a and t are the limits of the operator. There are many definitions for fractional derivatives. Three most frequently used ones are the Grunwald-Letnikov definition, the Riemann-Liouville, and the Caputo definitions.The Grunwald-Letnikov definition (GL) derivative with fractional-orderq is described by
(2)DtqaGLf(t)=lim⁡h→01hq∑j=0[(t-a)/h](-1)j(qj)f(t-jh),
where the symbol [•] means the integer part.The Riemann-Liouville (RL) definition of fractional derivatives is given by(3)DtqaRLf(t)=dndtn1Γ(n-q)∫atf(τ)(t-τ)q-n+1dτ,sfsdfsdfsdfdsfddjn-1<q<n,
where Γ(•) is the gamma function.The Caputo (C) fractional derivative is defined as follows:
(4)DtqaCf(t)=1Γ(n-q)∫at(t-τ)n-q-1f(n)(τ)dτ,fdsfsdfsdfdsdfddddln-1<q<n.It is well known that the initial conditions for the fractional differential equations with Caputo derivatives take the same form as for the integer-order ones, which is very suitable for practical problems [25]. Therefore, we will use the Caputo definition for the fractional derivatives in this paper.


### 2.2. Numerical Algorithms

Obtaining numerical solutions of a fractional differential equation is not easily compared with that for ordinary differential equations. There are two approximation methods which can frequently be used to numerical computations on chaos and bifurcations with fractional differential equations. One is an improved version of Adams-Bashforth-Moulton algorithm based on the predictor-correctors scheme [26], which is a time-domain approach. The other is a method, known as frequency domain approximation [27], based on numerical analysis of fractional-order systems in the frequency domain.Simulations for fractional-order systems using the time domain methods are complicated and due to long memory characteristics of these systems require a very long simulation time but on the other hand, it is more accurate [28]. Therefore, we employ the improved predictor-corrector algorithm for fractional-order differential equations in this paper.In order to get the approximate solution of a fractional-order chaotic system by the improved predictor-corrector algorithm, the following equation is considered:(5)dqxdtq=f(t,x),0≤t≤Tx(k)(0)=x0(k),k=0,1,…,⌈q⌉-1,
where ⌈q⌉ is just the value q rounded up to the nearest integer and x(k) is the ordinary kth derivative of x. Formula (5) is equivalent to the Volterra integral equation
(6)x(t)=∑k=0⌈q⌉-1x0(k)tkk!+1Γ(q)∫0t(t-τ)q-1f(τ,x(τ))dτ.
Now, for the sake of simplicity, we assume that we are working on a uniform grid {tn=nh:n=0,1,…,M} with some integer M and set h=T/M. Using the standard quadrature techniques for the integral in (6), denote g(τ)=f(τ,x(τ)); the integral is replaced by the trapezoidal quadrature formula at point tn+1(7)∫0tn+1(tn+1-τ)q-1g(τ)dτ≈∫0tn+1(tn+1-τ)q-1g~n+1(τ)dτ,
where g~n+1 is the piecewise linear interpolant for g with nodes and knots chosen at the tj, (j=0,1,…,n+1). After some elementary calculations, the right hand side of (7) gives
(8)∫0tn+1(tn+1-τ)q-1g~n+1(τ)dτ=hqq(q+1)∑j=0n+1αj,n+1g(tj).
And if we use the product rectangle rule, the right hand of (7) can be written as
(9)∫0tn+1(tn+1-τ)q-1g~n+1(τ)dτ=∑j=0n+1βj,n+1g(tj),
where
(10)αj,n+1={nq+1-(n-q)(n+1)q,fdsfddddddddj=0(n-j+2)q+1+(n-j)q+1-2(n-j+1)q+1,fdfdfdfdfdfdfdfdfdfdfdfdffff1≤j≤n,1,fsdfdfdsfdffdddfdfdfdfdffddj=n+1,βj,n+1=hqq((n-j+1)q-(n-j)q),1≤j≤n.
Then the predictor and corrector formulae for solving (6) are given, respectively, by
(11)xhP(tn+1)=∑k=0⌈q⌉-1x0(k)tn+1kk!+1Γ(q)∑j=0nβj,n+1f(tj,xh(tj)),xh(tn+1)=∑k=0⌈q⌉-1x0(k)tn+1kk!+hqΓ(q+2)f(tn+1,xhP(tn+1))jkjkjkkj+hqΓ(q+2)∑j=0nαj,n+1f(tj,xh(tj)).
The approximation accuracy of scheme (11) is O(hmin⁡{2,1+q}).


## 2.1. Definition

Fractional calculus is a generalization of integration and differentiation to a noninteger-order integrodifferential operatorDtqa which is defined by
(1)Dtqa={dqdtqR(q)>01R(q)=0∫at(dτ)-qR(q)<0,
where q is the derivative order which can be a complex number and R(q) is the real part of q. The numbers a and t are the limits of the operator. There are many definitions for fractional derivatives. Three most frequently used ones are the Grunwald-Letnikov definition, the Riemann-Liouville, and the Caputo definitions.The Grunwald-Letnikov definition (GL) derivative with fractional-orderq is described by
(2)DtqaGLf(t)=lim⁡h→01hq∑j=0[(t-a)/h](-1)j(qj)f(t-jh),
where the symbol [•] means the integer part.The Riemann-Liouville (RL) definition of fractional derivatives is given by(3)DtqaRLf(t)=dndtn1Γ(n-q)∫atf(τ)(t-τ)q-n+1dτ,sfsdfsdfsdfdsfddjn-1<q<n,
where Γ(•) is the gamma function.The Caputo (C) fractional derivative is defined as follows:
(4)DtqaCf(t)=1Γ(n-q)∫at(t-τ)n-q-1f(n)(τ)dτ,fdsfsdfsdfdsdfddddln-1<q<n.It is well known that the initial conditions for the fractional differential equations with Caputo derivatives take the same form as for the integer-order ones, which is very suitable for practical problems [25]. Therefore, we will use the Caputo definition for the fractional derivatives in this paper.


## 2.2. Numerical Algorithms

Obtaining numerical solutions of a fractional differential equation is not easily compared with that for ordinary differential equations. There are two approximation methods which can frequently be used to numerical computations on chaos and bifurcations with fractional differential equations. One is an improved version of Adams-Bashforth-Moulton algorithm based on the predictor-correctors scheme [26], which is a time-domain approach. The other is a method, known as frequency domain approximation [27], based on numerical analysis of fractional-order systems in the frequency domain.Simulations for fractional-order systems using the time domain methods are complicated and due to long memory characteristics of these systems require a very long simulation time but on the other hand, it is more accurate [28]. Therefore, we employ the improved predictor-corrector algorithm for fractional-order differential equations in this paper.In order to get the approximate solution of a fractional-order chaotic system by the improved predictor-corrector algorithm, the following equation is considered:(5)dqxdtq=f(t,x),0≤t≤Tx(k)(0)=x0(k),k=0,1,…,⌈q⌉-1,
where ⌈q⌉ is just the value q rounded up to the nearest integer and x(k) is the ordinary kth derivative of x. Formula (5) is equivalent to the Volterra integral equation
(6)x(t)=∑k=0⌈q⌉-1x0(k)tkk!+1Γ(q)∫0t(t-τ)q-1f(τ,x(τ))dτ.
Now, for the sake of simplicity, we assume that we are working on a uniform grid {tn=nh:n=0,1,…,M} with some integer M and set h=T/M. Using the standard quadrature techniques for the integral in (6), denote g(τ)=f(τ,x(τ)); the integral is replaced by the trapezoidal quadrature formula at point tn+1(7)∫0tn+1(tn+1-τ)q-1g(τ)dτ≈∫0tn+1(tn+1-τ)q-1g~n+1(τ)dτ,
where g~n+1 is the piecewise linear interpolant for g with nodes and knots chosen at the tj, (j=0,1,…,n+1). After some elementary calculations, the right hand side of (7) gives
(8)∫0tn+1(tn+1-τ)q-1g~n+1(τ)dτ=hqq(q+1)∑j=0n+1αj,n+1g(tj).
And if we use the product rectangle rule, the right hand of (7) can be written as
(9)∫0tn+1(tn+1-τ)q-1g~n+1(τ)dτ=∑j=0n+1βj,n+1g(tj),
where
(10)αj,n+1={nq+1-(n-q)(n+1)q,fdsfddddddddj=0(n-j+2)q+1+(n-j)q+1-2(n-j+1)q+1,fdfdfdfdfdfdfdfdfdfdfdfdffff1≤j≤n,1,fsdfdfdsfdffdddfdfdfdfdffddj=n+1,βj,n+1=hqq((n-j+1)q-(n-j)q),1≤j≤n.
Then the predictor and corrector formulae for solving (6) are given, respectively, by
(11)xhP(tn+1)=∑k=0⌈q⌉-1x0(k)tn+1kk!+1Γ(q)∑j=0nβj,n+1f(tj,xh(tj)),xh(tn+1)=∑k=0⌈q⌉-1x0(k)tn+1kk!+hqΓ(q+2)f(tn+1,xhP(tn+1))jkjkjkkj+hqΓ(q+2)∑j=0nαj,n+1f(tj,xh(tj)).
The approximation accuracy of scheme (11) is O(hmin⁡{2,1+q}).


## 3. The Chebyshev Polynomials

In general, random parameters in some engineering structures are bounded in nature. So the probability density function (PDF for short) model for the bounded random variables is taken to be an arch-like PDF of the form [29]:
(12)ρ(u)={(2π)1-u2|u|≤10|u|>1.Figure1 shows the arch-like PDF curve for bounded random variable u. As the orthogonal polynomial basis for this kind of PDF, the second kind of Chebyshev polynomials is the only choice. The general expansion for the second kind of Chebyshev polynomials can be described as follows:
(13)Hn(u)=∑k=0[n/2](-1)k(n-k)!k!(n-2k)!(2u)n-2k.
Then, we have
(14)H0(u)=1,H1(u)=2u,H2(u)=4u2-1,H3(u)=8u3-4u,H4(u)=16u4-12u2+1.⋯.
The recurrent formula for the Chebyshev polynomials of the second kind is
(15)uHn(u)=12[Hn-1(u)+Hn+1(u)].
The orthogonality of Chebyshev polynomial of the second kind can be expressed as
(16)∫-112π1-u2Hi(u)Hj(u)du={1(i=j)0(i≠j).
Equation (16) represents a weighted orthogonal relationship. The weighting function is just the same as that for the PDF of the random parameter which obeys arch distribution, ρ(u) in (12), and the left-hand side of (16) can be regarded as the expectation of the product Hi(u)Hj(u). Owing to the orthogonality of Chebyshev polynomial, any measurable function f(u)∈L2 can be expanded as
(17)f(u)=∑i=0∞ciHi(u),
where ci=∫-11ρ(u)f(u)Hi(u)du. This expansion is called an orthogonal decomposition of random function f(u), which is the theoretical base of orthogonal decomposition methods.Figure 1
The arch-like PDF curve for random variableu.


## 4. The Stochastic Fractional-Order VDP System

In this paper, we consider the fractional-order VDP system with a random parameter. The system can be described by the following differential equations:(18)Dqx=yDqy=μy-ax2y+bx,
where q is the fractional order and x, y are state variables. a, b are deterministic parameters and μ is a random parameter, and it can be expressed as
(19)μ=μ-+δu,
where μ- and δ/2 are the mean value and standard deviation of μ, respectively, and δ is regarded as the intensity of μ. u is a random variable defined on [-1,1].It is well known that the stochastic function space, which is constituted with the responses of nonlinear dynamic systems with random parameters, has proved to be a Hilbert space with respect to convergence in the mean square. Therefore, the fractional-order system (18) with random parameters can be reduced into the equivalent deterministic system by using of the orthogonal polynomial expansion. It follows from the orthogonal polynomial approximation that the responses of system (18) can be expressed approximately by the following series:
(20)x(t,u)=∑i=0Nxi(t)Hi(u)y(t,u)=∑i=0Nyi(t)Hi(u),
where Hi(u) represents the ith Chebyshev polynomial and N is the largest order of the polynomials. Then substituting (20) into the system (18), we can get
(21)Dq(∑i=0Nxi(t)Hi(u))=-(∑i=0Nyi(t)Hi(u))Dq(∑i=0Nyi(t)Hi(u))=μ-(∑i=0Nyi(t)Hi(u))+δu(∑i=0Nyi(t)Hi(u))-a(∑i=0Nxi(t)Hi(u))2(∑i=0Nyi(t)Hi(u))+b(∑i=0Nxi(t)Hi(u)),
and with the aid of the recurrent formulas of Chebyshev polynomial, all the nonlinear terms in the system can be further reduced into linear combination of Li(u) as follows:
(22)(∑i=0Nxi(t)Hi(u))2(∑i=0Nyi(t)Hi(u))2=∑i=03NSi(t)Hi(u),
where Si(t) can be obtained with the help of the computer algebraic systems, such as Maple or MATLAB. Meanwhile, the term δu(∑i=0Nyi(t)Hi(u)) can be reduced to
(23)uδ(∑i=0Nyi(t)Hi(u))=δ(∑i=0Nyi(t)Hi-1(u)+Hi+1(u)2).
Substituting (22) and (23) into (21), then we can have
(24)Dq(∑i=0Nxi(t)Hi(u))=-(∑i=0Nyi(t)Hi(u))Dq(∑i=0Nyi(t)Hi(u))=μ-(∑i=0Nyi(t)Hi(u))+δ(∑i=0Nyi(t)Hi-1(u)+Hi+1(u)2)-a∑i=0NSi(t)Hi(u)+b(∑i=0Nxi(t)Hi(u)).
Multiply both sides of (24) by Hi(u), i=0,1,2,…,N in sequence and take expectation with respect to u. Then we can finally obtain the equivalent deterministic system. Remember that if N→∞ in (24), then its left is strictly equivalent to the response of stochastic system (18). Otherwise, (24) is just approximately valid with a minimal residual error. Due to the requirement of computational precision, we take N=4 in the following numerical analysis. Therefore, we get the equivalent deterministic system approximately as
(25)Dqx0=-y0Dqy0=-aS0+bx0+μ-y0+12δy1Dqx1=-y1Dqy1=-aS1+bx1+μ-y1+12δ(y0+y2)Dqx2=-y2Dqy2=-aS2+bx2+μ-y2+12δ(y1+y3)Dqx3=-y3Dqy3=-aS3+bx3+μ-y3+12δ(y2+y4)Dqx4=-y4Dqy4=-aS4+bx4+μ-y4+12δy3.
We can get the numerical solutions xi, yi, (i=0,1,2,3,4) of system (25) by the improved predictor-corrector algorithm. Therefore, the approximate random responses of the stochastic system (18) can be expressed as
(26)x(t,u)≈∑i=04xi(t)Hi(u)y(t,u)≈∑i=04yi(t)Hi(u).
The ensemble mean responses (EMRs for short) of the stochastic fractional-order system are described as
(27)E[x(t,u)]≈∑i=04xi(t)E[Hi(u)]=x0(t)E[y(t,u)]≈∑i=04yi(t)E[Hi(u)]=y0(t).
Especially for u=0, namely, μ=μ-, this sample system is significant for reference. It is usually called a mean parameter system, of which the responses (SRMs for short) may be approximated as
(28)x(t,0)≈∑i=04xi(t)Hi(0)=x0(t)-x2(t)+x4(t)y(t,0)≈∑i=04yi(t)Hi(0)=y0(t)-y2(t)+y4(t).


## 5. Hopf Bifurcations

A Hopf bifurcation occurs when a periodic solution or limit cycle, surrounding an equilibrium point, arises or goes away as a parameter varies. In this section, the Hopf bifurcation for the stochastic fractional-order VDP system will be studied with theoretical and numerical computations. Obviously, both the stochastic VDP and equivalent deterministic systems (18) and (25) only have one equilibrium point, namely, the origin (0,0). The Jacobian matrix J for the linearized system evaluated at the origin is
(29)J=[0-100000000bμ-012δ000000000-100000012δ0bμ-012δ000000000-1000000012δbμ-012δ000000000-1000000012δbμ-012δ000000000-1000000012δbμ-].
The corresponding characteristic equation is
(30)f(λ)=a0λ10+a1λ9+a2λ8+a3λ7+a4λ6+a5λ5+a6λ4+a7λ3+a8λ2+a9λ+a10,
where a0,a1,…,a10 are the coefficients of (30) and can be obtained by computing via Maple or MATLAB. According to the bifurcation theory of nonlinear dynamical systems and taking w as a bifurcation parameter [22], we can obtain the existence conditions of Hopf bifurcation for the system. (1)
The eigenvalues of the linearized system about the fixed point areα(w)±iβ(w).(2)
The real and image parts of the eigenvalues satisfyα(wc)=0, β(wc)>0, α′(wc)≠0, when w=wc.(3)
The real parts of other eigenvalues are not equal to zeros.For the high order of equivalent system (25), so it is hard to get the threshold value of bifurcation parameter directly. Therefore, we can get the threshold via the following theorem.Theorem 1.
Letf(λ)=a0λn+a1λn-1+⋯+an=0 be the characteristic equation evaluated at the equilibrium point, Δn the n-dimensional Routh-Hurwitz determinant, and ai, (i=0,1,…,n) the coefficient of the characteristic equation. If there exists a threshold of bifurcation parameter wc which conforms to the following terms, then the first two existence conditions of Hopf bifurcation are satisfied:(1)
Δn-1(wc)=0, Δn-2(wc)≠0, Δn-3(wc)≠0;(2)
ai(wc)>0, (0,1,…,n);(3)
Δn-1′(wc)≠0.According to Theorem1, the following parameter relations which can make the equality Δ9=0 are derived:
(31)μ=0,μ=±δ2,μ=±δ4,μ=±3δ2,μ=±3δ4,μ=δ4±3δ4,μ=-δ4±3δ4.
Then, (31) is substituted into the 8-dimensional and 7-dimensional Routh-Hurwitz determinants Δ8, Δ7, and only the following parameters can make Δ8≠0, Δ7≠0 satisfied:
(32)μ=±δ2μ=±3δ2.
Meanwhile, the relations which can make all the coefficients satisfied ai>0 (i=0,1,…,10) are
(33)μ=±3δ2.
We substitute the relations (33) into Δ9′(μ) and find that Δ9′(μ) is not equal to 0; therefore, according to Theorem 1, the first two existence conditions of Hopf bifurcation are satisfied. When (33) is substituted into the characteristic equation, all the eigenvalues of it can be obtained:
(34)λ1,2=±bi,λ3,4=δ4-3δ4±144δ2-23δ2-16bλ5,6=-3δ4±143δ2-16bλ7,8=-δ4-3δ4±144δ2-23δ2-16bλ9,10=-3δ2±123δ2-4b.
Obviously, when the intensity of random parameter δ≠0, the real parts of all the other eigenvalues are less than 0 except λ1 and λ2. The characteristic equation does not have zero solution for b≠0. Therefore, all the existence conditions of Hopf bifurcation are satisfied. When μc=-3δ/2, the Hopf bifurcation occurs around the equilibrium (0,0).The parameters are taken asa=1, b=1, and δ=0.1, and the derivative order is q=0.995. And initial conditions are [x0(0),…,x4(0)]=(0.2,0.2,0.2,0.2)T and [y0(0),…,y4(0)]=(-0.2,-0.2,-0.2,-0.2)T. According to the above analysis, we can get that the critical parameter value of the bifurcation is μc=-3δ/2≈-0.0866, which implies that if μ<μc, the system will converge to the equilibrium; otherwise, it will converge to a limit cycle. The parameter values are chosen as μ=-0.09 and u=-0.01, respectively. Then we obtain the time trajectories and phase diagrams for the EMR and SRM of the stochastic system by simulation computations; see Figures 2 and 3. It is clearly seen that the system converges to the fixed point (0,0) for μ=-0.09 and to a limit cycle for u=-0.01, which implies that the Hopf bifurcation occurs in the stochastic system (18). Meanwhile, the SRM and EMR of phase plots coincide well enough. It means that the Chebyshev polynomial approximation method is effective for the stochastic fractional-order system (18).The time trajectories and phase diagrams foru=-0.09.
(a)
The time trajectory of EMR
(b)
The phase diagram of EMR inx0-y0 plane
(c)
The time trajectory of SRM
(d)
The phase diagram of SRMThe time trajectories and phase diagrams foru=-0.01.
(a)
The time trajectory of EMR
(b)
The phase diagram of EMR inx0-y0 plane
(c)
The time trajectory of SRM
(d)
The phase diagram of SRMWhen the others parameters are fixed andμ=-0.02, we also found that Hopf bifurcation occurs with the variation of the derivative order. The stochastic system converges to the equilibrium when q≤0.98 and to a limit cycle for q≥0.99, which means that the Hopf bifurcation occurs when the derivative order q=0.99; see Figure 4, from which it is obtained that the derivative order can also be taken as a bifurcation parameter for fractional-order system.The phase diagrams for different values of orderq.
(a)
q
=
0.9
8
(b)
q
=
0.99The simulation time is taken ast=800, one integration step of time is taken as Δt=0.001, derivative order q=0.97, and divide the time interval [0,800] into 8 equal parts. The time trajectory and phase diagram for t=800 are shown in Figure 5, from which it appears that the x0 converges to a fixed point as t≥300. However, the time trajectory is zoomed in when t belongs to the interval [300,400]; see Figure 6(a), from which we can see that the trajectory gradually decreases but does not converge to a point directly. Secondly, the time trajectory is zoomed in when t∈[400,500]; see Figure 6(b). Meanwhile, we also obtain the trajectories for t∈[500,600], t∈[600,700], and t∈[700,800] which are shown in Figures 6(c)–6(e), respectively. Comparing these figures, it is seen that the amplitude of x0 decreases as the time t increases. Based on these observations, we can suggest that the EMRs of the stochastic system (18) have an extremely long transient process while converging to the equilibrium as t→∞.The time trajectory and phase diagram of EMR whent∈[0,800]. (a) The time trajectory. (b) The phase diagram in x0-y0 plane.
(a)
(b)The time trajectories of EMRx0 for different time intervals.
(a)
t
∈
[
300,400
]
(b)
t
∈
[
400,500
]
(c)
t
∈
[
500,600
]
(d)
t
∈
[
600,700
]
(e)
t
∈
[
700,800
]


## 6. Conclusions

In this paper, the Hopf bifurcation of a fractional-order VDP system with a random parameter is analyzed. An equivalent deterministic system is obtained with the approximation principle of Chebyshev orthogonal polynomials. Based on bifurcation theory of nonlinear dynamic systems, the critical parameter value of the bifurcation is obtained by theoretical analysis. Meanwhile, numerical simulations are performed to verify the theoretical results. Besides, the derivative order can also be taken as a bifurcation parameter, and the Hopf bifurcation is observed with the variation of the order.It is worth mentioning that the orthogonal polynomial approximation is an effective method for integer-order stochastic systems, but the universality of the method for fractional-order ones still remains an open question. Due to the absences of analytical methods, the threshold value of bifurcation when the derivative order is taken as a bifurcation parameter cannot be obtained by theoretical analysis but only the numerical experiments.



---

*Source: 835482-2014-06-09.xml*