A Fractional Model for the Single Stokes Pulse from the Nonlinear Optics

Mihai Ivan

Teacher Training Department, Geometry and Topology Seminary, West University of Timișoara, Timișoara, Romania.

DOI: 10.4236/jamp.2022.1010191                

Abstract

In this paper we refer to equations of motion for the single Stokes pulse from the nonlinear optics, called the Stokes pulse system. A fractional-order model with Caputo derivative associated to Stokes pulse system (called the fractional Stokes pulse system) is proposed. The existence and uniqueness of solution of initial value problem for this fractional system are proved. The dynamic behavior for a special fractional Stokes pulse system is investigated, including: the fractional stability, the stabilization problem using suitable linear controls and the numerical integration based on fractional Euler method.

Keywords

Fractional Stokes Pulse System, Fractional Stability, Fractional Euler Method, Numerical Integration

1. Introduction

The theory of fractional differential equations (i.e. fractional calculus) and its applications are based on non-integer order of derivatives and integrals [1] [2].

The use of fractional models has received a great degree of interest in a series of works due to its applications in different fields of science and engineering. For example, these models played an important role in applied mathematics [3], mathematical physics [4], theoretical and applied physics [5], study of biological systems [6], control processing [7], chaos synchronization [8] [9] and so on. The dynamics of fractional-order systems associated to dynamical systems (in particular, Hamilton-Poisson systems) have been studied by many researchers in the recent decades [10] [11]. Another series of works deals with the study of dynamical behaviors of classical and fractional differential systems on Lie groups, Lie algebroids and Leibniz algebroids [12] [13] [14].

In this paper we consider the single Stokes pulse system [15]. It is described by the following differential equations on $R^3$ :

$\{\begin{array}{l}{\stackrel{\dot{}}{u}}^1\left(t\right)=\left(l_2-l_3\right)u^2\left(t\right)u^3\left(t\right)+b_2{u}^3\left(t\right)-b_3{u}^2\left(t\right)\\ {\stackrel{\dot{}}{u}}^2\left(t\right)=\left(l_3-l_1\right)u^3\left(t\right)u^1\left(t\right)+b_3{u}^1\left(t\right)-b_1{u}^3\left(t\right)\\ {\stackrel{\dot{}}{u}}^3\left(t\right)=\left(l_1-l_2\right)u^1\left(t\right)u^2\left(t\right)+b_1{u}^2\left(t\right)-b_2{u}^1\left(t\right)\end{array}$ (1.1)

where $u^1,u^2,u^3$ are state variables, ${\stackrel{\dot{}}{u}}^i\left(t\right)=\text{d}{u}^i\left(t\right)/\text{d}t$, $l_i,b_i\in R$ for $i=\overline{1,3}$ are parameters and t is the time.

The Hamilton-Poisson system (1.1) has been studied from mechanical geometry point of view [16]. It is associated to this system, the general fractional Stokes pulse system. The aim of our paper is focused on the study of a certain type of the fractional Stokes pulse system.

This paper is structured as follows. The Stokes pulse system (2.6) is described in Section 2. In Section 3 we define the fractional Stokes pulse system (3.1). The existence and uniqueness of solutions of initial value problem for the fractional model (3.1) are discussed. Also, are proposed four types of fractional Stokes pulse systems which are physically inequivalent. From the four types of fractional Stokes pulse systems, we choose a subcase of the second type, called the special fractional Stokes pulse system (3.3). The Section 4 is dedicated to analyzing of asymptotic stability of equilibrium states for the fractional model (3.3). For stabilization problem of the system (3.3), we associate the fractional Stokes pulse system with controls, denoted by (4.2). In Propositions (4.3) - (4.6) are established sufficient conditions on parameters k and k₁ to control the chaos in the fractional system (4.2). Using the fractional Euler’s method, the numerical integration of the system (4.2) is presented in Section 5.

2. The Single Stokes Pulse as Hamilton-Poisson System

For details on Hamiltonian dynamics, see e.g. [17] [18] [19].

The equations of motion for the Stokes polarization parameters of a single optical beam propagating as a traveling wave in a nonlinear medium (the single Stokes pulse) from the nonlinear optics are described by using the Stokes vector u and a Hamiltonian function H.

The Stokes vector is defined using the Pauli spin matrices and is called the polarization parameters of the single Stokes pulse. The vector $u\in R^3$ is assumed to be expressed in a linear polarization basis [15]. Let be A the transition matrix from the canonical basis of $R^3$ to the polarization basis. Since A is symmetric, then it one can always transform to a polarization basis in which A has the diagonal form $W=diag\left(l_1,l_2,l_3\right)$, where $l_i,i=\overline{1,3}$ are the eigenvalues of A.

The Hamiltonian function H is determined by the Stokes vector u, the diagonal matrix W, and the constant vectors $a=\left(a_1,a_2,a_3\right)$ and $c=\left(c_1,c_2,c_3\right)$. Using the vectors a and c, we define the vector $b=\left(b_1,b_2,b_3\right)$ by:

$b=a+r\cdot c,r=\left|c\right|.$ (2.1)

The matrix W describes the self-induced ellipse rotation. The vectors a and c describes the effects of linear and nonlinear anisotropy, respectively.

In terms of the Stokes parameters (the components of the vector u), the Hamiltonian function $H\in C^{\infty }\left(R^3,R\right)$, $u\to H\left(u\right)$, is defined by [15]:

$H\left(u\right)=u\cdot \left(b+\frac{1}{2}W\cdot u\right),$ (2.2)

where $b={\left(b_1,b_2,b_3\right)}^{\text{T}}$ and $u=u^{\text{T}}$.

The diagonal matrix W and the choice of the vectors a and c generates the dynamics of the Stokes vector u with the frequence b.

In the coordinate system $O{u}^1{u}^2{u}^3$, the Hamiltonian function H defined by (2.2), is written as:

$H\left(u\right)=\frac{1}{2}\left[l_1{\left(u^1\right)}^2+l_2{\left(u^2\right)}^2+l_3{\left(u^3\right)}^2\right]+b_1{u}^1+b_2{u}^2+b_3{u}^3.$ (2.3)

The dynamics of a single Stokes pulse is written as Hamilton-Poisson system. More precisely, the Stokes pulse system is defined on the Lie-Poisson manifold $so{\left(3\right)}^{\ast }$ the dual of Lie algebra $SO\left(3\right)$ with the following bracket:

$\left\{f,g\right\}\left(u\right)=u\cdot \left(\frac{\partial f}{\partial u}\times \frac{\partial g}{\partial u}\right)\text{,}(\forall )f,g\in C^{\infty }\left(R^{\text{3}}\text{,}R\right)$ (2.4)

and the Hamiltonian function $H:so{\left(3\right)}^{\ast }\cong R^3\to R$ given by (2.3).

The dynamical system defined on $so{\left(3\right)}^{\ast }$ with Poisson bracket {., .} given by (2.4), enabling the equations of motion to be expressed in Hamiltonian form:

$\stackrel{\dot{}}{u}=\left\{u,H\right\},$ (2.5)

where $u\in R^3$, t is the time and H is the Hamiltonian function [20].

We determine the equations ${\stackrel{\dot{}}{u}}^i=\left\{u^i,H\right\},i=\overline{1,3}$ of the system (2.5). We have:

$\begin{array}{c}\left\{u^1,H\right\}=u\cdot \left(\frac{\partial u^1}{\partial u}\times \frac{\partial H}{\partial u}\right)=\left|\begin{array}{ccc}{u}^1& u^2& u^3\\ 1& 0& 0\\ l_1{u}^1+b_1& l_2{u}^2+b_2& l_3{u}^3+b_3\end{array}\right|\\ =\left(l_2-l_3\right)u^2{u}^3+b_2{u}^3-b_3{u}^2.\end{array}$

Then, the first equation of the system (2.5) is ${\stackrel{\dot{}}{u}}^1=\left(l_2-l_3\right)u^2{u}^3+b_2{u}^3-b_3{u}^2$.

Finally one obtains the following differential system on $R^3$ :

$\{\begin{array}{l}{\stackrel{\dot{}}{u}}^1=\left(l_2-l_3\right)u^2{u}^3+b_2{u}^3-b_3{u}^2\\ {\stackrel{\dot{}}{u}}^2=\left(l_3-l_1\right)u^3{u}^1+b_3{u}^1-b_1{u}^3\\ {\stackrel{\dot{}}{u}}^3=\left(l_1-l_2\right)u^1{u}^2+b_1{u}^2-b_2{u}^1\end{array}$ (2.6)

where the parameters $l_i,b_i\in R$ for $i=\overline{1,3}$ are connected with the nature of the material and the medium.

The system (2.6) is called theStokes pulse dynamical system.

Acording to [20] there are six types of Equation (2.6) which are physically inequivalent and which correspond to different types of optical media [16].

Proposition 2.1. The functions H given by (2.4) and $C\in C^{\infty }\left(R^3,R\right)$, defined by:

$u\to C\left(u\right)=\frac{1}{2}\left[{\left(u^1\right)}^2+{\left(u^2\right)}^2+{\left(u^3\right)}^2\right],$

are constants of the motion (first integrals)for the dynamics (2.6).

Proof. Indeed, we have

$\text{d}H/\text{d}t=l_1{u}^1{\stackrel{\dot{}}{u}}^1+l_2{u}^2{\stackrel{\dot{}}{u}}^2+l_3{u}^3{\stackrel{\dot{}}{u}}^3+b_1{\stackrel{\dot{}}{u}}^1+b_2{\stackrel{\dot{}}{u}}^2+b_3{\stackrel{\dot{}}{u}}^3=0.$

Also, we have $\begin{array}{c}\text{d}C/\text{d}t=u^1{\stackrel{\dot{}}{u}}^1+u^2{\stackrel{\dot{}}{u}}^2+u^3{\stackrel{\dot{}}{u}}^3\\ =u^1\left[\left(l_2-l_3\right)u^2{u}^3+b_2{u}^3-b_3{u}^2\right]+u^2\left[\left(l_3-l_1\right)u^3{u}^1+b_3{u}^1-b_1{u}^3\right]\\ +u^3\left[\left(l_1-l_2\right)u^1{u}^2+b_1{u}^2-b_2{u}^1\right]\\ =0.\end{array}$

Remark 2.1. By Proposition 2.1, it follows thatthe trajectories of motion of Stokes pulse dynamical system (2.6) are intersections of the surfaces:

$H=\text{constant}$ and $C=\text{constant}$.

3. The Fractional Stokes Pulse System

For basic knowledge on fractional calculus, one may refer to [21] [22].

In this paper we consider the fractional derivative operator $D_t^q$ with $q\in \left(0,1\right)$ to be Caputo’s derivative. This fractional derivative operator is often used in concrete applications.

Let $f\in C^{\infty }\left(R\right)$ and $q\in R,q>0$. The q-order Caputo differential operator [21], is described by $D_t^{q}f\left(t\right)=I^{m-q}{f}^{\left(m\right)}\left(t\right),q>0$, where $f^{\left(m\right)}\left(t\right)$ represents the m-order derivative of the function f, $m\in N^{\ast }$ is an integer such that

$m-1\le q\le m$ and $I^q$ is $I^{q}f\left(t\right)=\frac{1}{\Gamma \left(q\right)}{\displaystyle \underset{0}{\overset{t}{\int }}{\left(t-s\right)}^{q-1}f\left(s\right)\text{d}s}$, $q>0$,

where $\Gamma$ is the Euler Gamma function. If $q=1$, then $D_t^{1}f\left(t\right)=\text{d}f/\text{d}t$.

The Hamilton-Poisson system (2.6) is modeled by the following fractional differential equations:

$\{\begin{array}{l}{D}_q^t{u}^1\left(t\right)=\left(l_2-l_3\right)u^2\left(t\right)u^3\left(t\right)+b_2{u}^3\left(t\right)-b_3{u}^2\left(t\right)\\ D_q^t{u}^2\left(t\right)=\left(l_3-l_1\right)u^3\left(t\right)u^1\left(t\right)+b_3{u}^1\left(t\right)-b_1{u}^3\left(t\right),q\in \left(0,1\right)\\ D_q^t{u}^3\left(t\right)=\left(l_1-l_2\right)u^1\left(t\right)u^2\left(t\right)+b_1{u}^2\left(t\right)-b_2{u}^1\left(t\right)\end{array}$ (3.1)

where $u^i$ are the Stokes polarization parameters and $l_i,b_i\in R$ for $i=\overline{1,3}$.

The system (3.1) is called the fractional Stokes pulse system associated to (2.6).

The initial value problem of fractional model (3.1) can be represented in the following matrix form:

$D_t^{q}u\left(t\right)=Au\left(t\right)+u^1\left(t\right)A_{1}u\left(t\right)+u^2\left(t\right)A_{2}u\left(t\right)+u^3\left(t\right)A_{3}u\left(t\right),u\left(0\right)=u_0\text{,}$ (3.2)

where $0<q<1$, $u\left(t\right)={\left(u^1\left(t\right),u^2\left(t\right),u^3\left(t\right)\right)}^{\text{T}}$, $t\in \left(0,\tau \right)$ and

$\begin{array}{l}A=\left(\begin{array}{ccc}0& -b_3& b_2\\ b_3& 0& -b_1\\ -b_2& b_1& 0\end{array}\right),A_1=\left(\begin{array}{ccc}0& 0& 0\\ 0& 0& l_3-l_1\\ 0& 0& 0\end{array}\right),\\ A_2=\left(\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ l_1-l_2& 0& 0\end{array}\right),A_3=\left(\begin{array}{ccc}0& l_2-l_3& 0\\ 0& 0& 0\\ 0& 0& 0\end{array}\right)\text{.}\end{array}$

Proposition 3.1.The initial value problem of the fractional Stokes pulse system (3.1)has a unique solution.

Proof. Let $f\left(u\left(t\right)\right)=Au\left(t\right)+u^1\left(t\right)A_{1}u\left(t\right)+u^2\left(t\right)A_{2}u\left(t\right)+u^3\left(t\right)A_{3}u\left(t\right)$. It is obviously continuous and bounded on $D=\left\{u\in R^3|u^i\in \left[u_0^i-\delta ,u_0^i+\delta \right],i=\overline{1,3}\right\}$ for any $\delta >0$. We have $f\left(u\left(t\right)\right)-f\left(u_1\left(t\right)\right)=A\left(u\left(t\right)-u_1\left(t\right)\right)+x\left(t\right)+y\left(t\right)+z\left(t\right)$, where $x\left(t\right)=u^1\left(t\right)A_{1}u\left(t\right)-u_1^1\left(t\right)A_1{u}_1\left(t\right)$, $y\left(t\right)=u^2\left(t\right)A_{2}u\left(t\right)-u_1^2\left(t\right)A_2{u}_1\left(t\right)$ and $z\left(t\right)=u^3\left(t\right)A_{3}u\left(t\right)-u_1^3\left(t\right)A_3{u}_1\left(t\right)$,

Then

(1) $\left|f\left(u\left(t\right)\right)-f\left(u_1\left(t\right)\right)\right|\le \Vert A\Vert \cdot \left|u\left(t\right)-u_1\left(t\right)\right|+\left|x\left(t\right)\right|+\left|y\left(t\right)\right|+\left|z\left(t\right)\right|$, where $\Vert \cdot \Vert$ and denote $\left|\cdot \right|$ matrix norm and vector norm, respectively.

It is easy to see that $x\left(t\right)=\left(u^1\left(t\right)-u_1^1\left(t\right)\right)A_{1}u\left(t\right)+u_1^1\left(t\right)A_1\left(u\left(t\right)-u_1\left(t\right)\right)$. Then

(2) $\left|x\left(t\right)\right|\le \Vert A_1\Vert \left(\left|u\left(t\right)\right|+\left|u_1^1\left(t\right)\right|\right)\left|u\left(t\right)-u_1\left(t\right)\right|$.

Similarly, we prove that

(3) $\left|y\left(t\right)\right|\le \Vert A_2\Vert \left(\left|u\left(t\right)\right|+\left|u_1^2\left(t\right)\right|\right)\left|u\left(t\right)-u_1\left(t\right)\right|$.

(4) $\left|z\left(t\right)\right|\le \Vert A_3\Vert \left(\left|u\left(t\right)\right|+\left|u_1^3\left(t\right)\right|\right)\left|u\left(t\right)-u_1\left(t\right)\right|$.

According to (2)-(4), the relation (1) becomes

(5) $\begin{array}{c}\left|f\left(u\left(t\right)\right)-f\left(u_1\left(t\right)\right)\right|\le (\Vert A\Vert +\Vert A_1\Vert \left(\left|u\left(t\right)\right|+\left|u_1^1\left(t\right)\right|\right)+\Vert A_2\Vert \left(\left|u\left(t\right)\right|+\left|u_1^2\left(t\right)\right|\right)\\ +\Vert A_3\Vert \left(\left|u\left(t\right)\right|+\left|u_1^3\left(t\right)\right|\right))\cdot \left|u\left(t\right)-u_1\left(t\right)\right|\end{array}$.

We have $\Vert A\Vert =\sqrt{2\left(b_1^2+b_2^2+b_3^2\right)}:=\theta$, $\Vert A_1\Vert =\left|l_3-l_1\right|$, $\Vert A_2\Vert =\left|l_1-l_2\right|$, $\Vert A_3\Vert =\left|l_2-l_3\right|$. If $\omega =\max \left\{\left|l_3-l_1\right|,\left|l_1-l_2\right|,\left|l_2-l_3\right|\right\}$, then $\Vert A_i\Vert \le \omega$ for $i=\overline{1,3}$. From the relation (5) we deduce that

(6) $\begin{array}{c}\left|f\left(u\left(t\right)\right)-f\left(u_1\left(t\right)\right)\right|\le (\theta +\omega \left(\left|u\left(t\right)\right|+\left|u_1^1\left(t\right)\right|\right)+\omega \left(\left|u\left(t\right)\right|+\left|u_1^2\left(t\right)\right|\right)\\ +\omega \left(\left|u\left(t\right)\right|+\left|u_1^3\left(t\right)\right|\right))\cdot \left|u\left(t\right)-u_1\left(t\right)\right|\end{array}$.

Given that inequalities $\left|u\left(t\right)\right|+\left|u_1^i\left(t\right)\right|\le 2\left|u_0\right|+\delta$, $i=\overline{1,3}$, are valid, the relation (6) becomes

(7) $\left|f\left(u\left(t\right)\right)-f\left(u_1\left(t\right)\right)\right|\le L\cdot \left|u\left(t\right)-u_1\left(t\right)\right|$,

where $L=\theta +3\omega \left(2\left|u_0\right|+\delta \right)>0$.

The inequality (7) shows that $f\left(u\left(t\right)\right)$ satisfies a Lipschitz condition. Using Theorems 1 and 2 in [23], it follows that the system (3.1) has a unique solution.

As with the nonlinear dynamics generated by the Stokes pulse system (2.6) there are six types of fractional Equation (3.1) which are physically inequivalent and which correspond to different types of optical media.

In this section we will refer to the types of fractional Stokes pulse systems for which the parameters $l_i\in R,i=\overline{1,3}$ meet the following condition $l_1\ne l_2\ne l_3\ne l_1$.

In this context there are the following four types of fractional Stokes pulse systems:

Type 1. $b=\left(0,0,0\right)$ and $l_1\ne l_2\ne l_3\ne l_1$ ;

Type 2. $b=\left(0,b_2,0\right)$, $b_2\ne 0$ and $l_1\ne l_2\ne l_3\ne l_1$ ;

Type 3. $b=\left(b_1,0,b_3\right)$, $b_1{b}_3\ne 0$ and $l_1\ne l_2\ne l_3\ne l_1$ ;

Type 4. $b=\left(b_1,b_2,b_3\right)$, $b_1{b}_2{b}_3\ne 0$ and $l_1\ne l_2\ne l_3\ne l_1$.

As an example, the fractional system corresponding to type 2, is given by:

$\{\begin{array}{l}{D}_q^t{u}^1\left(t\right)=\left(l_2-l_3\right)u^2\left(t\right)u^3\left(t\right)+b_2{u}^3\left(t\right)\\ D_q^t{u}^2\left(t\right)=\left(l_3-l_1\right)u^3\left(t\right)u^1\left(t\right),q\in \left(0,1\right)\\ D_q^t{u}^3\left(t\right)=\left(l_1-l_2\right)u^1\left(t\right)u^2\left(t\right)-b_2{u}^1\left(t\right)\end{array}$ (3.3)

where $b_2\in R^{\ast }$ and $l_i\in R,i=\overline{1,3}$ such that $l_1\ne l_2\ne l_3\ne l_1$.

The system (3.3) is called the special fractional Stokes pulse system. It is determined by a single nonzero component of vector b and contains the three nonlinear terms of the system (3.1) (given in general form).

For the system (3.3) we introduce the following notations:

$f_1\left(u\right)=\left(l_2-l_3\right)u^2{u}^3+b_2{u}^3,f_2\left(u\right)=\left(l_3-l_1\right)u^3{u}^1,f_3\left(u\right)=\left(l_1-l_2\right)u^1{u}^2-b_2{u}^1.$ (3.4)

Proposition 3.2. The equilibrium states of the special fractional Stokes pulse system (3.3) are given as the union of the following three families:

$\begin{array}{l}{E}_1:=\left\{e_1^m=\left(0,m,0\right)\in R^3|m\in R\right\},\\ E_2:=\left\{e_2^m=\left(m,\frac{b_2}{l_1-l_2},0\right)\in R^3|m\in R\right\},\\ E_3:=\left\{e_3^m=\left(0,\frac{b_2}{l_3-l_2},m\right)\in R^3|m\in R\right\}.\end{array}$

Proof. The equilibrium states are solutions of the equations $f_i\left(u\right)=0,i=\overline{1,3}$, where $f_i,i=\overline{1,3}$ are given by (3.4).

Remark 3.1. If in the fractional model (3.3) we take $q=1$, then one obtains the system for integer-order derivative which corresponds to type 4 $ of the dynamics (1.1). For this dynamical system, the nonlinear stability and the problem of existence of periodic solutions are studied, see Theorems 2.4, 3.5-3.7 [16].

4. Asymptotic Stability of the Special Fractional Stokes Pulse System (3.3)

Let us we present the study of asymptotic stability of equilibria for the fractional system (3.3). Finally, we will discuss how to stabilize the unstable equilibrium states of the system (3.3) via fractional order derivative. For this study we apply the Matignon’s test [7].

In the follows we will use the notations:

$\alpha :=l_2-l_3,\beta :=l_3-l_1,\gamma :=l_1-l_2.$ (4.1)

With the notations (4.1), the Jacobian matrix associated to system (3.3) is:

$J\left(u\right)=\left(\begin{array}{ccc}0& \alpha u^3& \alpha u^2+b_2\\ \beta u^3& 0& \beta u^1\\ \gamma u^2-b_2& \gamma u^1& 0\end{array}\right).$

Proposition 4.1. ( [7]) Let $u_e$ be an equilibrium state of system (3.3) and $J\left(u_e\right)$ be the Jacobian matrix $J\left(u\right)$ evaluated at $u_e$.

(i) $u_e$ is locally asymptotically stable, iff all eigenvalues of the matrix $J\left(u_e\right)$ satisfy:

$\left|\arg \left(\lambda \left(J\left(u_e\right)\right)\right)\right|>\frac{q\pi }{2}.$

(ii) $u_e$ is locally stable, iff either it is asymptotically stable, or the critical eigenvalues of $J\left(u_e\right)$ which satisfy $\left|\arg \left(\lambda \left(J\left(u_e\right)\right)\right)\right|=\frac{q\pi }{2}$ have geometric multiplicity one.

Proposition 4.2. The equilibrium states $e_i^m,i=\overline{1,3}$ are unstable $(\forall )q\in \left(0,1\right)$ .

Proof. The characteristic polynomial of the matrix $J\left(e_1^m\right)=\left(\begin{array}{ccc}0& 0& m\alpha +b_2\\ 0& 0& 0\\ m\gamma -b_2& 0& 0\end{array}\right)$ is $p_1\left(\lambda \right)=\det \left(J\left(e_1^m\right)-\lambda I\right)=-\lambda \left[{\lambda }^2-\left(\alpha \gamma m^2-b_2\left(\gamma -\alpha \right)m-b_2^2\right)\right]$ . For $m=0$, the characteristic polynomials of the matrix $J\left(e_0\right)$ is $p_0\left(\lambda \right)=-\lambda \left({\lambda }^2+b_2^2\right)$.

The characteristic polynomials of matrices $J\left(e_2^m\right)$ and $J\left(e_3^m\right)$ are the following:

$p_2\left(\lambda \right)=-\lambda \left({\lambda }^2-\beta \gamma m^2\right)$ and $p_3\left(\lambda \right)=-\lambda \left({\lambda }^2-\alpha \beta m^2\right)$ .

The equations $p_0\left(\lambda \right)=0$ and $p_i\left(\lambda \right)=0,i=\overline{1,3}$ have the root ${\lambda }_1=0$. Since $\arg \left({\lambda }_1\right)=0<\frac{q\pi }{2}$ for all $q\in \left(0,1\right)$, by Proposition 4.1 follows that the equilibrium states $e_0$ and $e_i^m,i=\overline{1,3}$ are unstable for all $q\in \left(0,1\right)$.

In the case when $u_e$ is a unstable equilibrium state of the fractional system (3.3), we associate to (3.3) a new fractional system, called the special fractional Stokes pulse system with (external)controls and given by:

$\{\begin{array}{l}{D}_q^t{u}^1=\alpha u^2{u}^3+b_2{u}^3+k_1{u}^1\\ D_q^t{u}^2=\beta u^3{u}^1+k{u}^2,q\in \left(0,1\right)\\ D_q^t{u}^3=\gamma u^1{u}^2-b_2{u}^1+k_1{u}^3\end{array}$ (4.2)

where $\alpha ,\beta ,\gamma$ are given in (4.1) and $k,k_1\in R$ are controls.

If one selects the parameters $k,k_1$ which then make the eigenvalues of the Jacobian matrix of fractional model (3.3) satisfy one of the conditions from Proposition 3.1, then its trajectories asymptotically approaches the unstable equilibrium state $u_e$ in the sense that ${\lim }_{t\to \infty }\Vert u\left(t\right)-u_e\Vert =0$, where $\Vert \cdot \Vert$ is the Euclidean norm.

The Jacobian matrix of the fractional model (4.2) with the controls $k,k_1$ is

$J\left(u,k,k_1\right)=\left(\begin{array}{ccc}{k}_1& \alpha u^3& \alpha u^2+b_2\\ \beta u^3& k& \beta u^1\\ \gamma u^2-b_2& \gamma u^1& k_1\end{array}\right).$

Proposition 4.3. Let be the fractional Stokes pulse system (4.2)with the controls $k,k_1\in R^{\ast }$.

(i) If $k<0,k_1<0$ , then $e_0$ is asymptotically stable $(\forall )q\in \left(0,1\right)$ ;

(ii) If $k<0,k_1>0$ and $q_0=\frac{2}{\pi }\arctan \frac{\left|b_2\right|}{k_1}$ , then:

(1) $e_0$ is asymptotically stable $(\forall )q\in \left(0,q_0\right)$ and it is stable for $q=q_0$ .

(2) $e_0$ is unstable $(\forall )q\in \left(q_0,1\right)$ .

(iii) If $k>0$ and $k_1\in R^{\ast }$ , then $e_0$ is unstable $(\forall )q\in \left(0,1\right)$ .

Proof. The characteristic polynomial of the Jacobian matrix $J\left(e_0,k,k_1\right)$ is $p_0\left(\lambda ,k,k_1\right)=-\left(\lambda -k\right)\left[{\left(\lambda -k_1\right)}^2+b_2^2\right]$. The roots of the equation $p_0\left(\lambda ,k,k_1\right)=0$ are ${\lambda }_1=k$, ${\lambda }_{2,3}=k_1\pm i{b}_2$.

(i) We suppose $k<0$ and $k_1<0$ . In this case we have $\mathrm{Re}\left({\lambda }_i\right)<0$ for $i=\overline{1,3}$. Since $\left|\arg \left({\lambda }_i\right)\right|=\pi >\frac{q\pi }{2},i=\overline{1,3}$ for all $q\in \left(0,1\right)$, by Proposition 4.1(i), it implies that $e_0$ is asymptotically stable for all $k,k_1\in R^{\ast }$.

(ii) We suppose $k<0$ and $k_1>0$ . In this case we have ${\lambda }_1<0$ and $\mathrm{Re}\left({\lambda }_{2,3}\right)>0$. Applying Proposition 4.1(i), $e_0$ is locally asymptotically stable, for $0<q<q_0$, where $q_0=\frac{2}{\pi }\arctan \frac{\left|b_2\right|}{k_1}$. If $q=q_0$, $e_0$ is stable. For $q_0<q<1$, $e_0$ is unstable $(\forall )k,k_1\in R^{\ast }$. Hence, the assertion (ii) holds.

(iii) We suppose $k>0$ and $k_1\in R^{\ast }$ . Since $J\left(e_0,k,k_1\right)$ has at least a positive eigenvalue, it follows that $e_0$ is unstable. Hence, (iii) holds, $(\forall )q\in \left(0,1\right)$.

Proposition 4.4. Let be the fractional Stokes pulse system (4.2)with the controls $k,k_1\in R^{\ast }$, $e_1^m=\left(0,m,0\right)$ and ${\Delta }_1=\left(\alpha m+b_2\right)\left(\gamma m-b_2\right)$ .

1. Let ${\Delta }_1<0$ and $q\in \left(0,1\right)$ .

(i) If $k<0$ and $k_1<0$ , then $e_1^m$ is asymptotically stable.

(ii) Let $k<0,k_1>0$ and $q_1=\frac{2}{\pi }\arctan \frac{\sqrt{-{\Delta }_1}}{k_1}$ .

(1) If $b_2>0,\frac{\alpha +\gamma }{\alpha \gamma }>0$ and $m\in \left(-\frac{b_2}{\alpha },\frac{b_2}{\gamma }\right)$ , then $e_1^m$ is asymptotically stable $(\forall )q\in \left(0,q_1\right)$ , stable for $q=q_1$ and unstable $(\forall )q\in \left(q_1,1\right)$ .

(2) If $b_2<0,\frac{\alpha +\gamma }{\alpha \gamma }>0$ and $m\in \left(\frac{b_2}{\gamma },-\frac{b_2}{\alpha }\right)$ , then $e_1^m$ is asymptotically stable $(\forall )q\in \left(0,q_1\right)$ , stable for $q=q_1$ and unstable $(\forall )q\in \left(q_1,1\right)$ .

(iii) If $k>0$ and $k_1\in R^{\ast }$ , then $e_1^m$ is unstable $(\forall )q\in \left(0,1\right)$ .

2. Let ${\Delta }_1>0$ and $q\in \left(0,1\right)$ .

(i) If $k<0,k_1<0$ and $k_1^2>{\Delta }_1$ .

(1) If $b_2>0,\frac{\alpha +\gamma }{\alpha \gamma }>0$ and $m\in \left(-\infty ,-\frac{b_2}{\alpha }\right)\cup \left(\frac{b_2}{\gamma },\infty \right)$ , then $e_1^m$ is asymptotically stable.

(2) If $b_2<0,\frac{\alpha +\gamma }{\alpha \gamma }>0$ and $m\in \left(-\infty ,\frac{b_2}{\gamma }\right)\cup \left(-\frac{b_2}{\alpha },\infty \right)$ , then $e_1^m$ is asymptotically stable.

(ii) Let $k<0,k_1<0,k_1^2\le {\Delta }_1$ and $\frac{\alpha +\gamma }{\alpha \gamma }>0$ . If $b_2>0$ , $m\in \left(-\infty ,-\frac{b_2}{\alpha }\right)\cup \left(\frac{b_2}{\gamma },\infty \right)$ or $b_2<0$ , $m\in \left(-\infty ,\frac{b_2}{\gamma }\right)\cup \left(-\frac{b_2}{\alpha },\infty \right)$, then $e_1^m$ is unstable.

Proof. The characteristic polynomial of (4.2) at $e_1^m$ is $J\left(e_1^m,k,k_1\right)=\left(\begin{array}{ccc}{k}_1& 0& \alpha m+b_2\\ 0& k& 0\\ \gamma m-b_2& 0& k_1\end{array}\right)$ whose characteristic polynomial is $p_1\left(\lambda ,k,k_1\right)=\det \left(J\left(e_1^m,k,k_1\right)-\lambda I\right)=-\left(\lambda -k\right)\left[{\left(\lambda -k_1\right)}^2-\left(\alpha m+b_2\right)\left(\gamma m-b_2\right)\right]$. The roots of the characteristic equation $p_1\left(\lambda ,k,k_1\right)=0$ are ${\lambda }_1=k$, ${\lambda }_{2,3}=k_1\pm \sqrt{{\Delta }_1}$, where ${\Delta }_1=\left(\alpha m+b_2\right)\left(\gamma m-b_2\right)$.

1. Case ${\Delta }_1<0$ and $q\in \left(0,1\right)$. We have the following two situations:

(1) if $b_2>0,\frac{\alpha +\gamma }{\alpha \gamma }>0$, then ${\Delta }_1<0$ for all $m\in \left(-\frac{b_2}{\alpha },\frac{b_2}{\gamma }\right)$ ;

(2) if $b_2<0,\frac{\alpha +\gamma }{\alpha \gamma }>0$, then ${\Delta }_1<0$ for all $m\in \left(\frac{b_2}{\gamma },-\frac{b_2}{\alpha }\right)$.

In this case, ${\lambda }_1=k$, ${\lambda }_{2,3}=k_1\pm i\sqrt{-{\Delta }_1}$.

(i) We suppose $k<0$ and $k_1<0$ . In this case we have ${\lambda }_1<0$ and $\mathrm{Re}\left({\lambda }_{2,3}\right)<0$. Since $\left|\arg \left({\lambda }_i\right)\right|=\pi >\frac{q\pi }{2},i=\overline{1,3}$ for all $q\in \left(0,1\right)$, by Proposition 4.1(i), it implies that $e_1^m$ is locally asymptotically stable for all $m\in R^{\ast }$.

(ii) (1)-(2). For $k<0$ and $k_1>0$, we have ${\lambda }_1<0$ and $\mathrm{Re}\left({\lambda }_{2,3}\right)>0$. Applying Proposition 4.1(i), $e_1^m$ is asymptotically stable, for $0<q<q_1$, where $q_1=\frac{2}{\pi }\arctan \frac{\sqrt{-{\Delta }_1}}{k_1}$. If $q=q_1$, then $e_1^m$ is stable. For $q_1<q<1$, $e_1^m$ is unstable.

(iii) Let $k>0$ and $k_1\in R^{\ast }$. Since ${\lambda }_1>0$, $J\left(e_1^m,k,k_1\right)$ has at least a positive eigenvalue and so $e_1^m$ is unstable. Hence, the assertions (i)-(iii) hold.

2. Case ${\Delta }_1>0$ and $q\in \left(0,1\right)$. Then ${\lambda }_1=k$, ${\lambda }_{2,3}=k_1\pm \sqrt{{\Delta }_1}$. We have the following two situations:

(1) if $b_2>0,\frac{\alpha +\gamma }{\alpha \gamma }>0$, then ${\Delta }_1>0$ for all $m\in \left(-\infty ,-\frac{b_2}{\alpha }\right)\cup \left(\frac{b_2}{\gamma },\infty \right)$ ;

(2) if $b_2<0,\frac{\alpha +\gamma }{\alpha \gamma }>0$, then ${\Delta }_1>0$ for all $m\in \left(-\infty ,\frac{b_2}{\gamma }\right)\cup \left(-\frac{b_2}{\alpha },\infty \right)$.

(i)-(ii) The eigenvalues ${\lambda }_i,i=\overline{1,3}$ are all negative if and only if $k<0,k_1<0$ and $k_1^2>{\Delta }_1$ . In these hypotheses it folows that $e_1^m$ is asymptotically stable. Also, if $k_1^2\le {\Delta }_1$, then $e_1^m$ is unstable. Therefore, the assertions (i)-(ii) hold.

Proposition 4.5. Let be the fractional Stokes pulse system (4.2)with the controls $k,k_1\in R^{\ast }$, $e_2^m=\left(m,\frac{b_2}{\gamma },0\right)$ and ${\Delta }_2={\left(k-k_1\right)}^2+4\beta \gamma m^2$ .

1. Let $\beta \gamma <0$ and $m\in R^{\ast }$ .

(i) Let ${\Delta }_2<0$ and $q\in \left(0,1\right)$ .

(1) If $k<0$ and $k_1<0$ , then $e_2^m$ is asymptotically stable $(\forall )m\in \left(-\infty ,\frac{k-k_1}{2\sqrt{-\beta \gamma }}\right)\cup \left(\frac{k_1-k}{2\sqrt{-\beta \gamma }},\infty \right)$ .

(2) If $k_1<0,k>k_1$ and $k+k_1<0$ , then $e_2^m$ is asymptotically stable $(\forall )m\in \left(-\infty ,\frac{k_1-k}{2\sqrt{-\beta \gamma }}\right)\cup \left(\frac{k-k_1}{2\sqrt{-\beta \gamma }},\infty \right)$ .

(3) If $k_1<0,k>0,k+k_1>0$ , $q_2=\frac{2}{\pi }\arctan \frac{\sqrt{-{\Delta }_2}}{k+k_1}$ and $m\in \left(-\infty ,\frac{k_1-k}{2\sqrt{-\beta \gamma }}\right)\cup \left(\frac{k-k_1}{2\sqrt{-\beta \gamma }},\infty \right)$ , then $e_2^m$ is asymptotically stable $(\forall )q\in \left(0,q_2\right)$ , stable for $q=q_2$ and unstable $(\forall )q\in \left(q_2,1\right)$ .

(ii) Let ${\Delta }_2\ge 0$ and $q\in \left(0,1\right)$ .

(1) If $k<k_1<0$ and $k_1<0$ , then $e_2^m$ is asymptotically stable $(\forall )m\in \left[\frac{k-k_1}{2\sqrt{-\beta \gamma }},\frac{k_1-k}{2\sqrt{-\beta \gamma }}\right]$ .

(2) If $k_1<0,k>k_1$ and $k+k_1<0$ , then $e_2^m$ is asymptotically stable $(\forall )m\in \left[\frac{k_1-k}{2\sqrt{-\beta \gamma }},\frac{k-k_1}{2\sqrt{-\beta \gamma }}\right]$ .

(3) If $k_1<0,k>0,k+k_1>0$ or $k_1>0$ , $k\in R^{\ast }$ , then $e_2^m$ is unstable.

2. Let $\beta \gamma >0$, $m\in R^{\ast }$ and $q\in \left(0,1\right)$ .

(i) Let $k<0,k_1<0$ and ${\left(k+k_1\right)}^2>{\Delta }_2>0$ .

(1) If $k<k_1$ , $m\in \left[\frac{k-k_1}{2\sqrt{-\beta \gamma }},\frac{k_1-k}{2\sqrt{-\beta \gamma }}\right]$ , then $e_2^m$ is asymptotically stable.

(2) If $k>k_1$ , $m\in \left[\frac{k_1-k}{2\sqrt{-\beta \gamma }},\frac{k-k_1}{2\sqrt{-\beta \gamma }}\right]$ , then $e_2^m$ is asymptotically stable.

(ii) Let $k_1<0,k+k_1>0,{\left(k+k_1\right)}^2\le {\Delta }_2$ or $k_1>0,k\in R^{\ast }$ . If $k<k_1$ , $m\in \left[\frac{k-k_1}{2\sqrt{-\beta \gamma }},\frac{k_1-k}{2\sqrt{-\beta \gamma }}\right]$ or $k>k_1,$ $m\in \left[\frac{k_1-k}{2\sqrt{-\beta \gamma }},\frac{k-k_1}{2\sqrt{-\beta \gamma }}\right]$ , then $e_2^m$ is unstable.

Proof. The Jacobian matrix of (4.2) at $e_2^m$ is $J\left(e_2^m,k,k_1\right)=\left(\begin{array}{ccc}{k}_1& 0& \frac{\alpha b_2}{\gamma }+b_2\\ 0& k& \beta m\\ 0& \gamma m& k_1\end{array}\right)$

whose characteristic polynomial is $p_2\left(\lambda ,k,k_1\right)=-\left(\lambda -k_1\right)\left[{\lambda }^2-\left(k+k_1\right)\lambda +k{k}_1-\beta \gamma m^2\right]$. The roots of the equation $p_2\left(\lambda ,k,k_1\right)=0$ are ${\lambda }_1=k_1$, ${\lambda }_{2,3}=\frac{\left(k+k_1\right)\pm \sqrt{{\Delta }_2}}{2}$, where ${\Delta }_2={\left(k-k_1\right)}^2+4\beta \gamma m^2$.

1. Case $\beta \gamma <0$ .

(i) Case ${\Delta }_2<0$ and $q\in \left(0,1\right)$. Then ${\lambda }_1=k_1$, ${\lambda }_{2,3}=\frac{\left(k+k_1\right)\pm i\sqrt{-{\Delta }_2}}{2}$. We have ${\Delta }_2<0$ if and only if $m\in \left(-\infty ,\frac{k-k_1}{2\sqrt{-\beta \gamma }}\right)\cup \left(\frac{k_1-k}{2\sqrt{-\beta \gamma }},\infty \right)$ when $k<k_1$ or $m\in \left(-\infty ,\frac{k_1-k}{2\sqrt{-\beta \gamma }}\right)\cup \left(\frac{k-k_1}{2\sqrt{-\beta \gamma }},\infty \right)$ when $k>k_1$.

(1) We suppose $k<k_1<0$. In this case we have ${\lambda }_1<0$ and $\mathrm{Re}\left({\lambda }_{2,3}\right)<0$. Since $\left|\arg \left({\lambda }_i\right)\right|=\pi >\frac{q\pi }{2},i=\overline{1,3}$ for all $q\in \left(0,1\right)$, by Proposition 4.1(i), it implies that $e_2^m$ is locally asymptotically stable for all $m\in \left(-\infty ,\frac{k-k_1}{2\sqrt{-\beta \gamma }}\right)\cup \left(\frac{k_1-k}{2\sqrt{-\beta \gamma }},\infty \right)$.

(2) We suppose $k_1<0,k>k_1$ and $k+k_1<0$. In this case we have ${\lambda }_1<0$ and $\mathrm{Re}\left({\lambda }_{2,3}\right)<0$. Applying the same reasoning as in the case (i)(1), one obtains that $e_2^m$ is locally asymptotically stable for all $m\in \left(-\infty ,\frac{k_1-k}{2\sqrt{-\beta \gamma }}\right)\cup \left(\frac{k-k_1}{2\sqrt{-\beta \gamma }},\infty \right)$.

(3) We suppose $k_1<0,k>0$ and $k+k_1>0$. In this case, ${\Delta }_2<0$ if and only if $m\in \left(-\infty ,\frac{k_1-k}{2\sqrt{-\beta \gamma }}\right)\cup \left(\frac{k-k_1}{2\sqrt{-\beta \gamma }},\infty \right)$. Then ${\lambda }_1<0$ and $\mathrm{Re}\left({\lambda }_{2,3}\right)>0$. Applying Proposition 4.1(i), $e_2^m$ is locally asymptotically stable, for $0<q<q_2$, where $q_2=\frac{2}{\pi }\arctan \frac{\sqrt{-{\Delta }_2}}{k+k_1}$. If $q=q_2$, then $e_2^m$ is stable. For $q_2<q<1$, $e_2^m$ is unstable.

(ii) Case ${\Delta }_2\ge 0$ and $q\in \left(0,1\right)$ . Then ${\lambda }_1=k_1$, ${\lambda }_{2,3}=\frac{\left(k+k_1\right)\pm \sqrt{{\Delta }_2}}{2}$. We have ${\Delta }_2\ge 0$ if and only if $m\in \left[\frac{k-k_1}{2\sqrt{-\beta \gamma }},\frac{k_1-k}{2\sqrt{-\beta \gamma }}\right]$ when $k<k_1$ or $m\in \left[\frac{k_1-k}{2\sqrt{-\beta \gamma }},\frac{k-k_1}{2\sqrt{-\beta \gamma }}\right]$ when $k>k_1$.

(1) and (2). In these cases, if $k<0$ and $k_1<0$, then ${\lambda }_2+{\lambda }_3<0$ and ${\lambda }_2{\lambda }_3=k{k}_1-\beta \gamma m^2>0$. It follows ${\lambda }_i<0,i=\overline{1,3}$. Then, $e_2^m$ is locally asymptotically stable.

(3) Let $k_1<0,k>0$ and $k+k_1>0$ or $k_1>0$ and $k\in R^{\ast }$ . Then $J\left(e_2^m,k,k_1\right)$ has at least a positive eigenvalue and so $e_2^m$ is unstable. Hence, the assertion (ii) holds.

2. Case $\beta \gamma >0$ , $m\in R^{\ast }$ and $q\in \left(0,1\right)$ . In this case ${\Delta }_2\ge 0$.

(i) We suppose $k_1<0$. We have ${\lambda }_2<0$ and ${\lambda }_3<0$ if and only if ${\lambda }_2+{\lambda }_3<0$ and ${\lambda }_2{\lambda }_3>0$. Then $k+k_1<0$ and ${\left(k+k_1\right)}^2-{\Delta }_2>0$. It follows ${\lambda }_i<0,i=\overline{1,3}$ for all $m\in R^{\ast }$ such that ${\left(k+k_1\right)}^2>{\Delta }_2$. Hence, $e_2^m$ is locally asymptotically stable.

(ii) We suppose $k_1<0,k+k_1>0$ or $k_1>0,k\in R^{\ast }$. Then, $J\left(e_2^m,k,k_1\right)$ has at least a positive eigenvalue and so $e_2^m$ is unstable. Therefore, the assertion (ii) holds.

Proposition 4.6. Let be the fractional Stokes pulse system (4.2)with the controls $k,k_1\in R^{\ast }$, $e_3^m=\left(0,-\frac{b_2}{\alpha },m\right)$ and ${\Delta }_3={\left(k-k_1\right)}^2+4\alpha \beta m^2$ .

1. Let $\alpha \beta <0$ and $m\in R^{\ast }$ .

(i) Let ${\Delta }_3<0$ and $q\in \left(0,1\right)$ .

(1) If $k<k_1<0$ , then $e_3^m$ is asymptotically stable $(\forall )m\in \left(-\infty ,\frac{k-k_1}{2\sqrt{-\alpha \beta }}\right)\cup \left(\frac{k_1-k}{2\sqrt{-\alpha \beta }},\infty \right)$ .

(2) If $k_1<0,k>k_1$ and $k+k_1<0$ , then $e_3^m$ is asymptotically stable $(\forall )m\in \left(-\infty ,\frac{k_1-k}{2\sqrt{-\alpha \beta }}\right)\cup \left(\frac{k-k_1}{2\sqrt{-\alpha \beta }},\infty \right)$ .

(3) If $k_1<0,k>0,k+k_1>0$ , $q_3=\frac{2}{\pi }\arctan \frac{\sqrt{-{\Delta }_3}}{k+k_1}$ and $m\in \left(-\infty ,\frac{k_1-k}{2\sqrt{-\alpha \beta }}\right)\cup \left(\frac{k-k_1}{2\sqrt{-\alpha \beta }},\infty \right)$ , then $e_2^m$ is asymptotically stable $(\forall )q\in \left(0,q_3\right)$ , stable for $q=q_3$ and unstable $(\forall )q\in \left(q_3,1\right)$ .

(ii) Let ${\Delta }_3\ge 0$ and $q\in \left(0,1\right)$ .

(1) If $k<k_1<0$ and $k_1<0$ , then $e_3^m$ is asymptotically stable $(\forall )m\in \left[\frac{k-k_1}{2\sqrt{-\alpha \beta }},\frac{k_1-k}{2\sqrt{-\alpha \beta }}\right]$ .

(2) If $k_1<0,k>k_1$ and $k+k_1<0$ , then $e_3^m$ is asymptotically stable $(\forall )m\in \left[\frac{k_1-k}{2\sqrt{-\alpha \beta }},\frac{k-k_1}{2\sqrt{-\alpha \beta }}\right]$ .

(3) If $k_1<0,k>0,k+k_1>0$ or $k_1>0$ , $k\in R^{\ast }$ , then $e_3^m$ is unstable.

2. Let $\alpha \beta >0$, $m\in R^{\ast }$ and $q\in \left(0,1\right)$ .

(i) Let $k<0,k_1<0$ and ${\left(k+k_1\right)}^2>{\Delta }_3>0$ .

(1) If $k<k_1$ , $m\in \left[\frac{k-k_1}{2\sqrt{-\alpha \beta }},\frac{k_1-k}{2\sqrt{-\alpha \beta }}\right]$ , then $e_3^m$ is asymptotically stable.

(2) If $k>k_1$ , $m\in \left[\frac{k_1-k}{2\sqrt{-\alpha \beta }},\frac{k-k_1}{2\sqrt{-\alpha \beta }}\right]$ , then $e_3^m$ is asymptotically stable.

(ii) Let $k_1<0,k+k_1>0,{\left(k+k_1\right)}^2\le {\Delta }_3$ or $k_1>0,k\in R^{\ast }$ . If $k<k_1$ , $m\in \left[\frac{k-k_1}{2\sqrt{-\alpha \beta }},\frac{k_1-k}{2\sqrt{-\alpha \beta }}\right]$ or $k>k_1$ , $m\in \left[\frac{k_1-k}{2\sqrt{-\alpha \beta }},\frac{k-k_1}{2\sqrt{-\alpha \beta }}\right]$, then $e_3^m$ is unstable.

Proof. The Jacobian matrix of (4.2) at $e_3^m$ is $J\left(e_3^m,k,k_1\right)=\left(\begin{array}{ccc}{k}_1& \alpha m& 0\\ \beta m& k& 0\\ -\frac{\gamma b_2}{\alpha }-b_2& 0& k_1\end{array}\right)$

Whose characteristic polynomial is $p_3\left(\lambda ,k,k_1\right)=-\left(\lambda -k_1\right)\left[{\lambda }^2-\left(k+k_1\right)\lambda +k{k}_1-\alpha \beta m^2\right]$. The roots of the equation $p_3\left(\lambda ,k,k_1\right)=0$ are ${\lambda }_1=k_1$, ${\lambda }_{2,3}=\frac{\left(k+k_1\right)\pm \sqrt{{\Delta }_3}}{2}$, where ${\Delta }_3={\left(k-k_1\right)}^2+4\alpha \beta m^2$.

1. Case $\alpha \beta <0$ .

(i) Case ${\Delta }_3<0$ and $q\in \left(0,1\right)$. Then ${\lambda }_1=k_1$, ${\lambda }_{2,3}=\frac{\left(k+k_1\right)\pm i\sqrt{-{\Delta }_3}}{2}$. We have ${\Delta }_3<0$ if and only if $m\in \left(-\infty ,\frac{k-k_1}{2\sqrt{-\alpha \beta }}\right)\cup \left(\frac{k_1-k}{2\sqrt{-\alpha \beta }},\infty \right)$ when $k<k_1$ or $m\in \left(-\infty ,\frac{k_1-k}{2\sqrt{-\alpha \beta }}\right)\cup \left(\frac{k-k_1}{2\sqrt{-\alpha \beta }},\infty \right)$ when $k>k_1$.

(1) We suppose $k<k_1<0$. In this case we have ${\lambda }_1<0$ and $\mathrm{Re}\left({\lambda }_{2,3}\right)<0$. Since $\left|\arg \left({\lambda }_i\right)\right|=\pi >\frac{q\pi }{2},i=\overline{1,3}$ for all $q\in \left(0,1\right)$, by Proposition 4.1(i), it implies that $e_3^m$ is locally asymptotically stable for all $m\in \left(-\infty ,\frac{k-k_1}{2\sqrt{-\alpha \beta }}\right)\cup \left(\frac{k_1-k}{2\sqrt{-\alpha \beta }},\infty \right)$.

(2) We suppose $k_1<0,k>k_1$ and $k+k_1<0$. In this case we have ${\lambda }_1<0$ and $\mathrm{Re}\left({\lambda }_{2,3}\right)<0$. Applying the same reasoning as in the case (i)(1), one obtains that $e_3^m$ is locally asymptotically stable for all $m\in \left(-\infty ,\frac{k_1-k}{2\sqrt{-\alpha \beta }}\right)\cup \left(\frac{k-k_1}{2\sqrt{-\alpha \beta }},\infty \right)$.

(3) We suppose $k_1<0,k>0$ and $k+k_1>0$. In this case, ${\Delta }_3<0$ if and only if $m\in \left(-\infty ,\frac{k_1-k}{2\sqrt{-\alpha \beta }}\right)\cup \left(\frac{k-k_1}{2\sqrt{-\alpha \beta }},\infty \right)$. Then ${\lambda }_1<0$ and $\mathrm{Re}\left({\lambda }_{2,3}\right)>0$. Applying Proposition 4.1(i), $e_3^m$ is locally asymptotically stable, for $0<q<q_3$, where $q_3=\frac{2}{\pi }\arctan \frac{\sqrt{-{\Delta }_3}}{k+k_1}$. If $q=q_3$, then $e_3^m$ is stable. For $q_3<q<1$, $e_3^m$ is unstable.

(ii) Case ${\Delta }_3\ge 0$ and $q\in \left(0,1\right)$ . Then ${\lambda }_1=k_1$, ${\lambda }_{2,3}=\frac{\left(k+k_1\right)\pm \sqrt{{\Delta }_3}}{2}$. We have ${\Delta }_3\ge 0$ if and only if $m\in \left[\frac{k-k_1}{2\sqrt{-\alpha \beta }},\frac{k_1-k}{2\sqrt{-\alpha \beta }}\right]$ when $k<k_1$ or $m\in \left[\frac{k_1-k}{2\sqrt{-\alpha \beta }},\frac{k-k_1}{2\sqrt{-\alpha \beta }}\right]$ when $k>k_1$.

(1) and (2). In these cases, if $k<0$ and $k_1<0$, then ${\lambda }_2+{\lambda }_3<0$ and ${\lambda }_2{\lambda }_3=k{k}_1-\alpha \beta m^2>0$. It follows ${\lambda }_i<0,i=\overline{1,3}$. Then, $e_3^m$ is locally asymptotically stable.

(3) Let $k_1<0,k>0$ and $k+k_1>0$ or $k_1>0$ and $k\in R^{\ast }$ . Then $J\left(e_3^m,k,k_1\right)$ has at least a positive eigenvalue and so $e_3^m$ is unstable. Hence, the assertion (ii) holds.

2. Case $\alpha \beta >0$ , $m\in R^{\ast }$ and $q\in \left(0,1\right)$ . In this case ${\Delta }_3\ge 0$.

(i) We suppose $k_1<0$. We have ${\lambda }_2<0$ and ${\lambda }_3<0$ if and only if ${\lambda }_2+{\lambda }_3<0$ and ${\lambda }_2{\lambda }_3>0$. Then $k+k_1<0$ and ${\left(k+k_1\right)}^2-{\Delta }_3>0$. It follows ${\lambda }_i<0,i=\overline{1,3}$ for all $m\in R^{\ast }$ such that ${\left(k+k_1\right)}^2>{\Delta }_3$. Hence, $e_3^m$ is locally asymptotically stable.

(ii) We suppose $k_1<0,k+k_1>0$ or $k_1>0,k\in R^{\ast }$. Then, $J\left(e_3^m,k,k_1\right)$ has at least a positive eigenvalue and so $e_3^m$ is unstable. Therefore, the assertion (ii) holds.

Example 4.1. (i) Let be the special fractional Stokes pulse system (4.2). We select $l_1=1,l_2=0.5,l_3=1.5$ and $b_2=1$. Then $\alpha =-1,\beta =\gamma =0.5$. We have $\frac{\alpha +\gamma }{\alpha \gamma }=1$ and ${\Delta }_1=\left(-m+1\right)\left(0.5m-1\right)$.

(i) Chosing $k=-0.15,k_1=-0.4$ and $m=1.6$, it follows that ${\Delta }_1=0.12$ and $k_1^2=0.16>\Delta {}_1$. According to Proposition 4.4, 2.(i)(1) it follows that the equilibrium state $e_1=\left(0,1.6,0\right)$ is asymptotically stable for $q=0.8$.

(ii) For $k=-1,k_1=-0.3$ and $m=1.6$, follows t ${\Delta }_1=0.12$ and $k_1^2=0.09<\Delta {}_1$.

The conditions of Proposition 4.4, 2.(i)(2) are achieved. Then $e_1=\left(0,1.6,0\right)$ is unstable for $q=0.75$.

Using Matlab, in Table 1 we give a set of values for the parameters $l_i,b_2,k,k_1,i=\overline{1,3}$, the equilibrium states and corresponding eigenvalues of special fractional Stokes pulse system (4.2).

Table 1. The controls $k,k_1$, equilibrium states $e_i^m$ and corresponding eigenvavues.

5. Numerical Integration of the Special Fractional Stokes Pulse System (4.2)

In this section we start with some mathematical preliminaries of the fractional Euler’s method for solving initial value problem for fractional differential equations.

Consider the following general form of the initial value problem (IVP) with Caputo derivative:

$D_t^{q}y\left(t\right)=f\left(t,y\left(t\right)\right),y\left(0\right)=y_0,t\in I=\left[0,T\right],T>0,$ (5.1)

where $y:I\to R^n,f:R^n\to R^n$ is a continuous nonlinear function and $q\in \left(0,1\right)$, represents the order of the derivative.

The right-hand side of the IVP (5.1) in considered examples is Lipschitz functions and the numerical method used in this works to integrate system (5.1) is the Fractional Euler’s method.

Since f is assumed to be continuous function, every solution of the initial value problem given by (5.1) is also a solution of the following Volterra fractional integral equation:

$y\left(t\right)=y\left(0\right)+I_t^{q}f\left(t,y\left(t\right)\right),$ (5.2)

where $I_t^q$ is the q-order Riemann-Liouville integral operator, which is expressed by:

$I_t^{q}f\left(t\right)=\frac{1}{\Gamma \left(q\right)}{\displaystyle \underset{0}{\overset{t}{\int }}{\left(t-s\right)}^{q-1}f\left(s,y\left(s\right)\right)\text{d}s},q>0,s\in \left[0,T\right].$ (5.3)

Moreover, every solution of (5.2) is a solution of the (IVP) (5.1).

To integrate the fractional Equation (5.1), means to find the solution of (5.2) over the interval $\left[0,T\right]$. In this context, a set of points $\left(t_j,y\left(t_j\right)\right)$ are produced which are used as approximated values. In order to achieve this approximation, the interval $\left[0,T\right]$ is partitioned into n subintervals $\left[t_j,t_{j+1}\right]$ each equal width

$h=\frac{T}{n},t_j=jh$ for $j=0,1,\cdots ,n$.

For the fractional-order q and $j=0,1,2,\cdots$, it computes an approximation denoted as $y_{j+1}$ for $y\left(t_{j+1}\right),j=0,1,2,\cdots$.

The general formula of the fractional Euler’s method for to compute the elements $y_j$, is:

$y_{j+1}=y_j+\frac{h^q}{\Gamma \left(q+1\right)}f\left(t_j,y\left(t_j\right)\right),t_{j+1}=t_j+h\text{,}j=0,1,\cdots ,n\text{.}$ (5.4)

For more details, see [24] [25].

For the numerical integration of the special fractional Stokes pulse system (4.2), we apply the fractional Euler method (FEM). For this, consider the following fractional differential equations:

$\{\begin{array}{l}{D}_t^q{u}^i\left(t\right)=F_i\left(u^1\left(t\right),u^2\left(t\right),u^3\left(t\right)\right),i=\overline{1,3},t\in \left(t_0,\tau \right],q\in \left(0,1\right)\\ u\left(t_0\right)=\left(u^1\left(t_0\right),u^2\left(t_0\right),u^3\left(t_0\right)\right),\end{array}$ (5.5)

where

$\{\begin{array}{l}{F}_1\left(u\left(t\right)\right)=\left(l_2-l_3\right)u^2\left(t\right)u^3\left(t\right)+b_2{u}^3\left(t\right)+k_1{u}^1\left(t\right)\\ F_2\left(u\left(t\right)\right)=\left(l_3-l_1\right)u^3\left(t\right)u^1\left(t\right)+k{u}^2\left(t\right),q\in \left(0,1\right)\\ F_3\left(u\left(t\right)\right)=\left(l_1-l_2\right)u^1\left(t\right)u^2\left(t\right)-b_2{u}^1\left(t\right)+k_1{u}^3\left(t\right)\end{array}$ (5.6)

where $b_2,k,k_1\in R^{\ast }$ and $l_i\in R,i=\overline{1,3}$ such that $l_1\ne l_2\ne l_3\ne l_1$.

Since the functions $F_i\left(u\left(t\right)\right),i=\overline{1,3}$ are continuous, then the initial value problem (5.5) is equivalent to system of Volterra integral equations, which is given as follows:

$u^i\left(t\right)=u^i\left(0\right)+I_t^q{F}_i\left(u^1\left(t\right),u^2\left(t\right),u^3\left(t\right)\right),i=\overline{1,3}\text{.}$ (5.7)

The system (5.7) is called the Volterra integral equations associated to special Stokes pulse system (4.2).

The problem for solving the system (5.5) is reduced to one of solving a sequence of systems of fractional equations in increasing dimension on successive intervals $\left[j,\left(j+1\right)\right]$.

For the numerical integration of the system (5.6) one can use the fractional Euler method (the formula (5.4), which is expressed as follows:

$u^i\left(j+1\right)=u^i\left(j\right)+\frac{h^q}{\Gamma \left(q+1\right)}{F}_i\left(u^1\left(j\right),u^2\left(j\right),u^3\left(j\right)\right),i=\overline{1,3}$ (5.8)

where $j=0,1,\cdots ,N$, $h=\frac{T}{N}$, $T>0$, $N>0$.

More precisely, the numerical integration of the fractional system (5.5) is given by:

$\{\begin{array}{l}{u}^1\left(j+1\right)=u^1\left(j\right)+h^q\frac{1}{\Gamma \left(q+1\right)}\left(\alpha u^2\left(j\right)u^3\left(j\right)+b_2{u}^3\left(j\right)+k_1{u}^1\left(j\right)\right)\\ u^2\left(j+1\right)=u^2\left(j\right)+h^q\frac{1}{\Gamma \left(q+1\right)}\left(\beta u^3\left(j\right)u^1\left(j\right)+k{u}^2\left(j\right)\right),q\in \left(0,1\right)\\ u^3\left(j+1\right)=u^3\left(j\right)+h^q\frac{1}{\Gamma \left(q+1\right)}\left(\gamma u^1\left(j\right)u^2\left(j\right)-b_2{u}^1\left(j\right)+k_1{u}^3\left(j\right)\right)\\ u^i\left(0\right)=u_e^i+\epsilon ,i=\overline{1,3},\end{array}$ (5.9)

where $\alpha :=l_2-l_3$, $\beta :=l_3-l_1$, $\gamma :=l_1-l_2$.

Using [21] [24], we have that the numerical algorithm given by (5.9) is convergent.

Example 5.1. Let us we present the numerical integration of the special fractional Stokes pulse system with controls which has considered in Example 4.1(i). For this we apply the algorithm (5.9) and software Maple. Then, in (5.9) we take: $l_1=1$, $l_2=0.5$, $l_3=1.5$, $b_2=1$, $k=-0.15$, and $k_1=-0.4$. It is known that the equilibrium state $e_1=\left(0,1.6,0\right)$ is asymptotically stable.

For the numerical simulation of solutions of the above fractional model we use the rutine Maple. spec-fract-Stokes-pulse-system-with-controls, denoted by [sp-fr.Stokes-pulse syst]. Applying this program for $h=0.01$, $\epsilon =0.01$, $u^1\left(0\right)=\epsilon$, $u^2\left(0\right)=1.6+\epsilon$, $u^3\left(0\right)=\epsilon$, $N=100$, $t=102$, one obtain the orbits $\left(n,u^1\left(n\right)\right)$, $\left(n,u^2\left(n\right)\right)$ $\left(n,u^3\left(n\right)\right)$ and $\left(u^1\left(n\right),u^2\left(n\right),u^3\left(n\right)\right)$, for $q=0.8$.

Finally, we present the rutine [sp-fr.Stokes pulse syst]:

# Fractional equations associated to Stokes pulse system for $q=0.8$

Du1/dt=(l2-l3)*u2*u3 + b2*u3 + k1* u1;

Du2/dt=(l3-l1)*u1*u3 + k* u2;

Du3/dt=(l1-l2)*u1*u2 - b2*u1 + k1* u3;

> with (plots):

> l1:=1.; l2:=0.5; l3:=1.5; alpha:=l2-l3; beta:=l3-l1; gamma:=l1-l2; b2:=1.; k:=-0.15; k1:= -0.4; q:=0.8; u1e:=0.; u2e:=1.6; u3e:=0.;

> with (stats):

> h:=0.01; epsilon:=0.01; n:=100:t:=n+2; u1:= array (0 .. n): u2:= array (0 .. n): u3:= array (0 .. n): u1[0]:=epsilon + u1e; u2[0]:=epsilon + u2e; u3[0]:=epsilon + u3e;

> for j from 1 by 1 to n do

> u1[j]:= u1[j-1] + h Ùq *(alpha* u2[j-1]*u3[j-1] + b2*u3[j-1] + k1* u1[j-1])/GAMMA(q+1);

u2[j]:= u2[j-1] + h Ùq *(beta* u1[j-1]*u3[j-1] + k* u2[j-1])/GAMMA(q+1);

u3[j]:= u3[j-1] + h Ùq *(gamma* u1[j-1]*u2[j-1] - b2*u1[j-1] + k1* u3[j-1])/GAMMA(q+1);

od:

> plot (seq([j,u1[j]], j = 0 .. n), style = point, symbol = point, scaling = UNCONSTRAINED);

plot (seq([j,u2[j]], j = 0 .. n), style = point, symbol = point, scaling = UNCONSTRAINED);

plot (seq([j,u3[j]], j = 0 .. n), style = point, symbol = point, scaling = UNCONSTRAINED);

pointplot 3d ( {seq([u1[j], [u2[j], [u3[j]], j = 0 .. n)}, style = point, symbol = point, scaling =

UNCONSTRAINED, color = red);

Remark 5.1. Appyling (5.9) and Maple for the numerical simulation of solutions of fractional model (4.2) for each set of values for parameters $l_1,l_2,l_3,b_2,k$ and $k_1$, given in the Table 1, it will be found that the results obtained are valid.

Conclusions. This paper presents the fractional Stokes pulse system (3.1) associated to system (2.6). The special fractional Stokes pulse system (3.3) was studied from fractional differential equations theory point of view: asymptotic stability, determining of sufficient conditions on parameters $k,k_1$ to control the chaos in the proposed fractional system and numerical integration of the fractional model (4.2). By choosing the right parameters $k$ and $k_1$ in the fractional model (4.2), this work offers a series of chaotic fractional differential systems. The other types of systems mentioned in the four types of fractional models will be studied in future works.

Acknowledgments

The author has very grateful to be reviewers for their comments and suggestions.

Conflicts of Interest

The author declares no conflicts of interest regarding the publication of this paper.

References

[1]	Podlubny, I. (1999) Fractional Differential Equations. Academic Press, Cambridge.
[2]	Kilbas, A.A., Srivastava, H.M. and Trujillo, J.J. (2006) Theory and Applications of Fractional Differential Equations. North-Holland Mathematics Studies, 204, 1-523.
[3]	Ahmed, H.M. (2018) Fractional Euler Method: An Effective Tool for Solving Fractional Differential Equations of Fractional Order. Journal of the Egyptian Mathematical Society, 26, 38-43. https://doi.org/10.21608/JOEMS.2018.9460
[4]	Agrawal, O.P. (2002) Solution for a Fractional Diffusion-Wave Equation Defined in a Bounded Domain. Nonlinear Dynamics, 29, 145-155. https://doi.org/10.1023/A:1016539022492
[5]	Hilfer, R. (2000) Applications of Fractional Calculus in Physics. World Scientific, Singapore. https://doi.org/10.1142/3779
[6]	Ahmed, E., El-Sayed, A.M.A. and El-Saka, H.A.A. (2007) Equilibrium Points, Stability and Numerical Solutions of Fractional-Order Predator-Prey and Rabies Models. Journal of Mathematical Analysis and Applications, 325, 542-553. https://doi.org/10.1016/j.jmaa.2006.01.087
[7]	Matignon, D. (1996) Stability Results for Fractional Differential Equations with Applications to Control Processing. Proceedings of the Computational Engineering in Systems and Applications, Vol. 2, Lille, July 1996, 963-968.
[8]	Zhou, T. and Li, C. (2005) Synchronization in Fractional-Order Differential Systems. Physica D Nonlinear Phenomena, 212, 111-125. https://doi.org/10.1016/j.physd.2005.09.012.
[9]	Bhalekar, S. (2012) Chaos Control and Synchronization in Fractional-Order Lorenz-Like System. International Journal of Differential Equations, 2012, Article ID: 623234, 16 p. https://doi.org/10.1155/2012/623234
[10]	Ahmed, E., El-Sayed, A.M.A. and El-Saka, H.A.A. (2006) On Some Routh-Hurwitz Conditions for Fractional Order Differential Equations and Their Applications in Lorenz, Rössler, Chua and Chen Systems. Physics Letters A, 358, 1-4. https://doi.org/10.1016/j.physleta.2006.04.087
[11]	Ivan, M. and Ivan, G. (2018) On the Fractional Euler Top System with Two Parameters. Journal of Modern Engineering Research, 8, 10-22.
[12]	Degeratu, M. and Ivan, M. (2006) Linear Connections on Lie Algebroids. Proceedings of the 5th Conference of Balkan Society of Geometers, Mangalia, 29 August-2 September 2005, 44-53.
[13]	Ivan, M., Ivan, G. and Opriș, D. (2009) Fractional Equations of the Rigid Body on the Pseudo-Orthogonal Group SO(2,1). International Journal of Geometric Methods in Modern Physics, 6, 1181-1192. https://doi.org/10.1142/S0219887809004168
[14]	Ivan, G., Ivan, M. and Opriș, D. (2007) Fractional Dynamical Systems on Fractional Leibniz Algebroids. Analele științifice ale Universității “Al. I. Cuza” din Iași (S.N.), West University of Timișoara, Timișoara, 222-234.
[15]	David, D., Holm, D.D., Tratnik, M.V. (1990) Hamiltonian Chaos in Nonlinear Optical Polarization Dynamics. Physics Reports, 187, 281-367. https://doi.org/10.1016/0370-1573(90)90063-8
[16]	Puta, M., Cașu, I. and Butur, M. (2002) Some Remarks on the Nonlinear Optical pol Arization Dynamics. Tensor. New Series, 63, 126-146.
[17]	Marsden, J.E. and Rațiu, T.S. (1999) Introduction to Mechanics and Symmetry. Springer, Berlin, 17.
[18]	Ivan, G. (2013) Geometrical and Dynamical Properties of General Euler Top System. Indian Journal of Pure and Applied Mathematics, 44, 77-93. https://doi.org/10.1007/s13226-013-0004-0
[19]	Pop, C., Aron, A., Galea, C., Ciobanu, M. and Ivan, M. (2009) Some Geometric Aspects in Theory of Lotka-Volterra System. Proceedings of the 11th WSEAS International Conference on Sustainability in Science Engineering, Iasi, 27 May 2009, 91-97.
[20]	Tratnik, M.V. and Sipe, J.E. (1987) Nonlinear Polarization Dynamics. I. The Single Pulse Equations. Physical Review A, 35, 2965-2975. https://doi.org/10.1103/PhysRevA.35.2965
[21]	Diethelm, K. (2010) The Analysis of Fractional Differential Equations: An Application-Oriented Exposition Using Differential Operators of Caputo Type. Springer, Berlin.
[22]	Băleanu, D., Diethelm, K., Scalas, E. and Trujillo J.J. (2012) Fractional Calculus Models and Numerical Methods. In: Luo, A.C.J., Ed., Series on Complexity, Nonlinearity and Chaos, Vol. 3, World Scientific, Singapore, 428. https://doi.org/10.1142/8180
[23]	Cheng, Z. (2011) Dynamic Analysis of a Fractional Order Rössler System. International Journal of Nonlinear Science, 12, 227-235.
[24]	Odibat, Z.M. and Momani, S. (2008) An Algorithm for the Numerical Solution of Differential Equations of Fractional Order. Journal of Applied Mathematics & Informatics, 26, 15-27. http://www.kcam.biz.
[25]	Danca, M.-F. (2017) Hidden Chaotic Attractors in Fractional-Order Systems. Nonlinear Dynamics, 89, 577-586. https://doi.org/10.1007/s11071-017-3472-7