Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order

*Kamal Shah, Eiman, Hammad Khalil, Rahmat Ali Khan and Thabet Abdeljawad*

#### **Abstract**

This paper is concerned about the study of shifted Jacobi polynomials. By means of these polynomials, we construct some operational matrices of fractional order integration and differentiations. Based on these matrices, we develop a numerical scheme for the boundary value problems of fractional order differential equations. The construction of the procedure is new one for the coupled systems of fractional order boundary value problems. In the proposed scheme, we obtain a simple but highly accurate systems of algebraic equations. These systems are easily soluble by means of Matlab or using Mathematica. We provide some examples to which the procedure is applied to verify the applicability of our proposed method.

**Keywords:** Jacobi polynomials, operational matrices, fractional order differential equations, coupled system, boundary conditions

#### **1. Introduction**

Recently, fractional order differential equations have gained much attention from the researchers of mathematics, physics, computer science and engineers. This is due to the large numbers of applications of fractional order differential equations in various fields of science and applied nature. The applications of fractional order differential equations are found in physics, mechanics, viscoelasticity, photography, biology, chemistry, fluid mechanics, image and signal processing phenomenons, etc.; for more details, see [1–4]. The researchers are giving more attention to fractional order differential equations, because in most cases the fractional order models are more accurate and reliable than classical order models. On the other hand, fractional order models have more degrees of freedom than integer order models. In most cases, the hereditary phenomena and memory description process are accurately modeled by means of fractional order derivatives and integrations as compare with integer order derivative and integrations. Due to these facts, researchers had given much attention to study the existence and uniqueness of positive solutions as well as multiplicity of solutions for fractional order differential equations and their systems. This area is very well explored, and large numbers of research articles can be found in literature; for details, see [5–10]. It is very difficult to obtain exact solutions for such type of differential equations in all cases always. However, the area involving numerical solutions is in initial stage. Therefore,

various numerical schemes have been developed in last few years by different researchers. The aim of these schemes was to achieve better accuracy. These schemes have their own merits and demerits. For example, in [11–13], the authors developed Adomain decomposition methods, Homotopy analysis methods, and Variational iterational methods to obtained good approximate solutions for certain fractional order differential equations. In recent years, spectral methods have got much attention as they were applied to solve some real-word problems of various fields of science and engineering. High accuracy was obtained for solving such problems [14]. Spectral method needs operational matrices for the numerical solutions, which have been constructed by using some polynomials, for example, in [15], the authors developed an operational matrix for shifted Legendre polynomials corresponding to fractional order derivative. In [16], the authors introduced operational matrix by using shifted Legendre polynomials corresponding to fractional order integrations. Similarly in [17, 18], authors constructed operational matrices for fractional order derivative by using Chebyshev and Jacobi polynomials. In all these cases, these matrices were applied to solve multiterms fractional order differential equations together with Tau-collocation method. In [19], Singh et al. derived an operational matrices for Berestein polynomials corresponding to fractional order integrations. Haar wavelet operational matrices were also developed, and some problems of fractional order differential equations were solved. In all these methods, collocation method was used together with these methods to obtain numerical solutions for fractional order differential equations. They only solved single problems.

Most of the physical and biological phenomena can be modeled in the form of coupled systems. For example, when two masses and two spring systems are modeled, it led to a coupled system of differential equations.

$$\begin{aligned} \mathbf{D}^{\alpha}X + \frac{\lambda}{m\_1} \mathbf{D}^{\nu\_1}X - \frac{\mu}{m\_1} \mathbf{D}^{\nu\_2}Y + \frac{\kappa\_1}{m\_1}X - \frac{\kappa\_2}{m\_1}Y &= f(x), \ 1 < a \le 2, \ 0 < \nu\_1, \nu\_2 \le 1, \\\mathbf{D}^{\beta}Y - \frac{\lambda}{m\_2} \mathbf{D}^{\alpha\_1}X + \frac{\lambda + \mu}{m\_2} \mathbf{D}^{\alpha\_2}Y - \frac{\kappa\_1}{m\_2}X + \frac{\kappa\_1 + \kappa\_2}{m\_2}Y &= g(x), \ 1 < \beta \le 2, \ 0 < \nu\_1, \nu\_2 \le 1, \end{aligned}$$

subject to the boundary conditions

$$\begin{aligned} X(\mathbf{0}) &= \mathbf{0}, & X(\mathbf{1}) &= \mathbf{0} \\ Y(\mathbf{0}) &= \mathbf{1}, & Y(\mathbf{1}) &= -\mathbf{1}. \end{aligned}$$

Many coupled systems can be found in physics and in biological model, chemistry, fluid dynamics, etc.; see [16, 20–27]. In fractional calculus, we generalize these systems by replacing the integer order by any fractional values that lie in certain range. All the above numerical methods were applied to solve singles problems, and to the best of our knowledge, very few articles are devoted to the numerical solutions of coupled systems or models. In all these cases, the authors have solved initial value problems of coupled systems of fractional order differential and partial differential equations; for details, see [15, 28] and the references therein. For example, in [15], authors have solved some initial value problems of coupled systems of fractional differential equations (FDEs) by using shifted Jacobi polynomials operational matrix method. However, the scheme is not properly applied yet for solving coupled systems of boundary value problems of fractional order differential equations. In this paper, we discuss the shifted Jacobi polynomial operational matrices methods to solve boundary value problems for a coupled systems of fractional order differential equations. In this scheme, we introduced a new matrix corresponding to boundary conditions, which is required for the approximate solutions. First, we solve single boundary value problems given by

*Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order DOI: http://dx.doi.org/10.5772/intechopen.102054*

$$\mathbf{D}^{\alpha}\mathbf{x}(t) + A\mathbf{D}^{\nu}\mathbf{x}(t) + B\mathbf{D}^{\alpha}\mathbf{x}(t) + \mathbf{C}\mathbf{x}(t) = f(t), t \in [0, T], A, B, \mathbf{C} \in \mathbb{R}, \tag{1}$$

subject the boundary conditions

$$\mathfrak{x}(\mathbf{0}) = \mathfrak{x}\_0, \quad \mathfrak{x}(T) = \mathfrak{x}\_1,$$

where *f t*ð Þ∈*C*½ � 0, *T* is a source term and 1< *α*≤2, 0 <*ν*, *ω*<1*:* Then, we extend our technique to solve coupled system of fractional differential equations (FDEs) given as

$$\begin{aligned} \mathbf{D}^{\alpha}\mathbf{x}(t) + A\_1 \mathbf{D}^{\nu\_1}\mathbf{x}(t) + B\_1 \mathbf{D}^{\nu\_2}\mathbf{y}(t) + C\_1 \mathbf{x}(t) + D\_1 \mathbf{y}(t) &= f(t), \ t \in [0, T], \\ \mathbf{D}^{\beta}\mathbf{y}(t) + A\_2 \mathbf{D}^{\alpha\_1}\mathbf{x}(t) + B\_2 \mathbf{D}^{\alpha\_2}\mathbf{y}(t) + C\_2 \mathbf{x}(t) + D\_2 \mathbf{y}(t) &= \mathbf{g}(t), \ t \in [0, T], \end{aligned} \tag{2}$$

subject to the boundary conditions given by

$$\begin{aligned} \varkappa(\mathbf{0}) &= \varkappa\_0, \ \varkappa(T) = \varkappa\_1, \\ \varkappa(\mathbf{0}) &= \mathcal{y}\_0, \ \jmath(T) = \mathcal{y}\_1, \end{aligned} \tag{3}$$

where 1<*α*, *β* ≤2, 0< *νi*,*ω<sup>i</sup>* <1ð Þ *i* ¼ 1, 2 *:* Further, *f t*ð Þ, *g t*ð Þ∈*C*½ � 0, *T* are the source term of the system (42), (43) and *x*0, *y*0, *x*1, *y*<sup>1</sup> are real constants, *Ai*, *Bi*,*Ci*, *Di*ð Þ *i* ¼ 1, 2 are any constants. We convert the differential equations/ systems of differential equations to a coupled systems of simple algebraic equations of Selvester or Laypunove type, which are easily soluble for unknown matrix required for approximations. In our proposed method, we do no need collation method, which is required for abovementioned methods. We applied our technique to some practical problems to verify the applicability of the proposed method.

Our article is organized as follows: Some fundamental concepts and results of fractional calculus and Jacobi polynomials needed are provided in sections 2. Section 3 contains some operational matrices for shifted Jacobi polynomials corresponding to fractional order derivative and integrations. Also in this section, we provide new operational matrix corresponding to the boundary conditions. In Section 4, we derive main procedure for the general class of boundary value problems of fractional order differential equations as well as for coupled system of FDEs. In Section 5, we apply the proposed method to some examples. Also we plot the images of these examples. In Section 6, we provide a short conclusion of the papers.

#### **2. Preliminaries**

In this section, we recall some basic definitions and results of fractional calculus, also we will give fundamental properties, definitions related to Jacobi polynomials as in [1, 4, 6–8, 14, 17, 18, 28–31].

**Definition 2.1.** The fractional integral of order *α* >0 of a function *z* : ð Þ! 0, ∞ *R* is defined by

$$\mathbf{I}^a z(t) = \frac{1}{\Gamma(a)} \int\_0^t \frac{z(s)}{\left(t - s\right)^{1 - a}} ds,\tag{4}$$

provided the integral converges at the right sides. Further, a simple and important property of **I** *<sup>α</sup>* is given by

$$\mathbf{I}^a \mathbf{z}^\delta = \frac{\Gamma(\delta + 1)}{\Gamma(\delta + a + 1)} \mathbf{z}^{a + \delta}. \tag{5}$$

**Definition 2.2.** *The Caputo fractional derivative of order α*>0 *of a function <sup>z</sup>*∈*Cn*½ � 0, 1 *is defined by*

$$\begin{aligned} \mathbf{D}\_{0+}^{a}z(t) &= \frac{\mathbf{1}}{\Gamma(n-a)} \int\_{0}^{t} (t-s)^{n-a-1} \frac{d^{n}}{dt^{n}} z(s)ds, & n-1 < a \le n, t > 0, \\\text{where } n = \lceil a \rceil + \mathbf{1}, \end{aligned} \tag{6}$$

provided that the right side is pointwise defined on 0, ð Þ ∞ . Also one of the important properties of fractional order derivative is given by

$$\mathbf{D}^a t^k = \frac{\Gamma(k+1)}{\Gamma(k-a+1)} t^{k-a}. \text{ Also for any constant } \mathbf{C}, \text{we have } \mathbf{D}^a \mathbf{C} = \mathbf{0}. \tag{7}$$

The following results are needed in the sequel. **Lemma 2.2.1.** [6], Let *α* >0 then

$$\mathbf{I}^a \mathbf{D}^a y(t) = y(t) + \mathbf{C}\_0 + \mathbf{C}\_1 t + \dots + \mathbf{C}\_{n-1} t^{a-n},\tag{8}$$

for arbitrary

$$\mathbf{C}\_i \in \mathbb{R}, \ i = 0, 1, 2, \dots, n - 1, \ n = [\alpha] + 1.$$

#### **2.1 The shifted Jacobi polynomials and its fundamental properties**

In this section, we provide basic properties of shifted Jacobi polynomials. The famous Jacobi polynomials *P*ð Þ *<sup>a</sup>*,*<sup>b</sup> <sup>i</sup>* ð Þ*z* are defined over the interval ½ � �1, 1 as given by

$$\begin{split} P\_i^{(a,b)}(\mathbf{z}) &= \frac{(a+b+2i-1)\left[a^2 - b^2 + z(a+b+2i-2)(a+b+2i-2)\right]}{2i(a+b+i)(a+b+2i-2)} P\_{i-1}^{(a,b)}(\mathbf{z}) \\ &- \frac{(a+i-1)(b+i-1)(a+b+2i)}{i(a+b+i)(a+b+2i-2)} P\_{i-2}^{(a,b)}(\mathbf{z}), \quad i=2,3,\quad \dots \, \end{split} \tag{6}$$
 
$$\text{where} \qquad P\_0^{(a,b)}(\mathbf{z}) = \mathbf{1}, \quad P\_1^{(a,b)}(\mathbf{z}) = \frac{a+b+2}{2}\mathbf{z} + \frac{a-b}{2}. \tag{9}$$

By means of the substitution *<sup>z</sup>* <sup>¼</sup> <sup>2</sup>*<sup>t</sup> <sup>ξ</sup>* � 1, we get a new version of the above polynomials defined in (9), which is called the shifted Jacobi polynomials over the interval 0, ½ � *ξ* . The analytical form of this shifted Jacobi polynomials is given by the relations

$$\begin{split} P^{(a,b)}\_{\xi,i}(t) &= \sum\_{k=0}^{i} \frac{(-1)^{i-k} \Gamma(i+b+1) \Gamma(i+k+a+b+1)}{\Gamma(k+b+1) \Gamma(i+a+b+1) \Gamma(i-k+1) \Gamma(k+1) \tilde{\xi}^{k}} t^{k}, \\ \text{where} \quad P^{(a,b)}\_{\xi,i}(0) &= (-1)^{i} \frac{\Gamma(i+b+1)}{\Gamma(b+1) \Gamma(i+1)}. \end{split} \tag{10}$$

The orthogonality conditions of the shifted Jacobi polynomials are given by

$$\begin{cases} \overset{\circ}{\rho}\_{\xi\dot{j}}^{(a,b)}(t)P\_{\xi,i}^{(a,b)}(t)W\_{\xi}^{(a,b)}(t)dt = \mathfrak{Q}\_{\xi\dot{j}}^{(a,b)}\delta\_{\vec{\mu}},\\ \text{where } \delta\_{\vec{\mu}} = \mathbf{1}, \quad \text{if } i = j, \quad \text{otherwise } \delta\_{\vec{\mu}} = \mathbf{0},\\ \text{and the weight function is given by } W\_{\xi}^{(a,b)}(t) = (\xi - t)^{a}t^{b}, \end{cases} \tag{11}$$

*Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order DOI: http://dx.doi.org/10.5772/intechopen.102054*

and

$$\mathfrak{Q}^{(a,b)}\_{\xi,j}(t) = \frac{\xi^{a+b+1} \Gamma(\lfloor j+a+1 \rfloor) \Gamma(\lfloor j+b+1 \rfloor)}{(2j+a+b+1) \Gamma(\lfloor j+1 \rfloor) \Gamma(\lfloor j+a+b+1 \rfloor)}. \tag{12}$$

Any function *f t*ð Þ such that *f* 2 ð Þ*t* is integrable over 0, ½ � *ξ* , can be approximated in terms of shifted Jacobi polynomials as given by

$$f(t) \simeq \sum\_{k=0}^{m} B\_j P\_{\xi,k}^{(a,b)}(t) = \mathbf{K}\_M^T \Phi\_M(t),\tag{13}$$

where the shifted Jacobi coefficient vector is denoted by K*<sup>M</sup>* and Φ*M*ð Þ*t* is *M* terms function vector. Also *M* ¼ *m* þ 1, when *m* ! ∞ the approximation converges to the exact value of the function. The coefficient *bk* can be calculated by using (11)–(13) as

$$B\_{j} = \frac{1}{\Omega\_{\xi\_{j}^{\circ}}^{(a,b)}} \int\_{0}^{\xi} W\_{\xi}^{(a,b)}(t) f(t) P\_{\xi,i}^{(a,b)}(t) dt, \; i = 0, 1, \ldots \tag{14}$$

#### **3. Operational matrices**

In this section, we provide some operational matrices for fractional order differentiation as well as for fractional order integration. Further, in this section, we give construction of new matrix for the boundary conditions.

**Theorem 3.1.** Let Φ*M*ð Þ*t* be the function vector, then the operational matrix of fractional order derivative *α* is given by

$$\mathbf{D}^{a}[\Phi\_{\mathcal{M}}(t)] \simeq \mathbf{H}^{(a)}\_{\mathcal{M}\times\mathcal{M}}\Phi\_{\mathcal{M}}(t),\tag{15}$$

where Hð Þ *<sup>α</sup> <sup>M</sup>*�*<sup>M</sup>* is the operational matrix of fractional order derivative *<sup>α</sup>* given by

$$\mathbf{H}\_{M\times M}^{(a)} = \begin{bmatrix} \mathbf{Y}\_{0,0,k} & \mathbf{Y}\_{0,1,k} & \cdots & \mathbf{Y}\_{0,j,k} & \cdots & \mathbf{Y}\_{0,m,k} \\ \mathbf{Y}\_{2,1,k} & \mathbf{Y}\_{2,2,k} & \cdots & \mathbf{Y}\_{2,r,k} & \cdots & \mathbf{Y}\_{1,m,k} \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \\ \mathbf{Y}\_{i,0,k} & \mathbf{Y}\_{i,1,k} & \cdots & \mathbf{Y}\_{ij,k} & \cdots & \mathbf{Y}\_{i,m,k} \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \\ \mathbf{Y}\_{m,0,k} & \mathbf{Y}\_{m,1,k} & \cdots & \mathbf{Y}\_{m,j,k} & \cdots & \mathbf{Y}\_{m,m,k} \end{bmatrix},\tag{16}$$

where

$$\begin{split} \Upsilon\_{i,j,k} &= \sum\_{k=0}^{i} \frac{(-1)^{i-k} \Gamma(i+b+1) \Gamma(i+k+a+b+1) \Gamma(k+1)}{\Gamma(k+b+1) \Gamma(i+a+b+1) \Gamma(i-k+1) \Gamma(k+1) \Gamma(k+1-a) \xi^{k}} \\ &\times \sum\_{l=0}^{j} \frac{(-1)^{j-l} (2j+a+b+1) \Gamma(j+1) \Gamma(j+l+a+b+1) \Gamma(k-alpha+l+b+1) \Gamma(a+1) \xi^{a}}{\Gamma(j+a+1) \Gamma(l+b+1) \Gamma(j-l+1) \Gamma(l+1) \Gamma(k-a+l+b+a+2)} \\ &\neq f' \quad \Upsilon\_{i,j,k} = 0, \quad for \ i < a. \end{split} \tag{17}$$

*Proof.* The proof of this theorem is same as given in [14, 17] in lemma 3ð Þ *:*2 and ð Þ <sup>3</sup>*:*<sup>4</sup> , respectively. □

**Theorem 3.2.** *Let* **Φ***M*ð Þ*t be the function vector corresponding to the shifted Jacobi polynomials and α*>0*, then the operational matrix corresponding to the fractional order integration is given by*

$$\mathbf{H}^{a}[\Phi\_{\mathbf{M}}(\mathbf{t})] \simeq \mathbf{H}^{(a)}\_{\mathbf{M}\times\mathbf{M}}\Phi\_{\mathbf{M}}(\mathbf{t}),\tag{18}$$

where Hð Þ *<sup>α</sup> <sup>M</sup>*�*<sup>M</sup>* is the operational matrix of fractional order integration *<sup>α</sup>* and is given by

$$\mathbf{H}\_{M\times M}^{(a)} = \begin{bmatrix} \mathbf{Y}\_{0,0,k} & \mathbf{Y}\_{0,1,k} & \cdots & \mathbf{Y}\_{0,j,k} & \cdots & \mathbf{Y}\_{0,m,k} \\ \mathbf{Y}\_{2,1,k} & \mathbf{Y}\_{2,2,k} & \cdots & \mathbf{Y}\_{2,r,k} & \cdots & \mathbf{Y}\_{1,m,k} \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \\ \mathbf{Y}\_{i,0,k} & \mathbf{Y}\_{i,1,k} & \cdots & \mathbf{Y}\_{ij,k} & \cdots & \mathbf{Y}\_{i,m,k} \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \\ \mathbf{Y}\_{m,0,k} & \mathbf{Y}\_{m,1,k} & \cdots & \mathbf{Y}\_{m,j,k} & \cdots & \mathbf{Y}\_{m,m,k} \end{bmatrix},\tag{19}$$

where

$$\begin{split} \Upsilon\_{i,j,k} &= \sum\_{k=0}^{i} \frac{(-1)^{i-k} \Gamma(i+b+1) \Gamma(i+k+a+b+1) \Gamma(k+1)}{\Gamma(k+b+1) \Gamma(i+a+b+1) \Gamma(i-k+1) \Gamma(k+1) \Gamma(k+1+a) \xi^{k}} \\ &\times \sum\_{l=0}^{j} \frac{(-1)^{j-l} (2j+a+b+1) \Gamma(j+1) \Gamma(j+l+a+b+1) \Gamma(k+a+l+b+1) \Gamma(a+1) \xi^{a}}{\Gamma(j+a+1) \Gamma(l+b+1) \Gamma(j-l+1) \Gamma(l+1) \Gamma(k+a+l+b+a+2)} . \end{split} \tag{20}$$

*Proof.* The proof of this theorem is available in [14]. Therefore, we omit the proof.

**Theorem 3.3** *Let* Φ*M*ð Þ*t be a function vector, and let ϕ*ð Þ*t be any function defined as <sup>ϕ</sup>*ðÞ¼ *<sup>t</sup> ct<sup>n</sup>*, *<sup>n</sup>* <sup>¼</sup> 0, 1, 2, … ,*c*<sup>∈</sup> *and u t*ðÞ¼ <sup>K</sup>*<sup>T</sup> <sup>M</sup>*Φ*M*ð Þ*t . Then*

$$
\phi(t)\_0 \mathbf{I}\_\xi^a u(t) = \mathbf{K}\_M^T \mathbf{Q}\_{M \times M}^{(c, \phi, a)} \Phi\_M(t), \tag{21}
$$

where Qð Þ *<sup>c</sup>*,*ϕ*,*<sup>α</sup> <sup>M</sup>*�*<sup>M</sup> is an operational matrix corresponding to some boundary value and is given by*

$$\mathbf{Q}\_{M\times M}^{(c,\phi,a)} = \begin{bmatrix} \mathbf{Y}\_{0,0,k} & \mathbf{Y}\_{0,1,k} & \cdots & \mathbf{Y}\_{0,j,k} & \cdots & \mathbf{Y}\_{0,m,k} \\ \mathbf{Y}\_{2,1,k} & \mathbf{Y}\_{2,2,k} & \cdots & \mathbf{Y}\_{2,r,k} & \cdots & \mathbf{Y}\_{1,m,k} \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \\ \mathbf{Y}\_{i,0,k} & \mathbf{Y}\_{i,1,k} & \cdots & \mathbf{Y}\_{ij,k} & \cdots & \mathbf{Y}\_{i,m,k} \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \\ \mathbf{Y}\_{m,0,k} & \mathbf{Y}\_{m,1,k} & \cdots & \mathbf{Y}\_{mj,k} & \cdots & \mathbf{Y}\_{m,m,k} \end{bmatrix},\tag{22}$$

where

$$\Upsilon\_{ij,k} = \sum\_{j=0}^{M} \Delta\_{i,k,a} B\_j P\_{\xi\_j^j}^{(a,b)}(t),\tag{23}$$

where

$$\Delta\_{i,k,a} = \sum\_{k=0}^{i} \frac{(-1)^{i-k} (i+b+1) \Gamma(i+a+b+1) c \xi^a}{\Gamma(k+b+1) \Gamma(i+a+b+1) \Gamma(i-k+1) \Gamma(k+a)} \tag{24}$$

*Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order DOI: http://dx.doi.org/10.5772/intechopen.102054*

and

$$B\_{j} = \sum\_{l=0}^{j} \frac{(-1)^{j-l} (2j+a+b+1) \Gamma(j+1) \Gamma(j+l+a+b+1) \Gamma(n+b+1) \Gamma(a+1) \xi^{n-l}}{\Gamma(j+a+1) \Gamma(l+b+1) \Gamma(j-l+1) \Gamma(l+1) \Gamma(n+a+b+1)}. \tag{25}$$

*Proof.* Consider the general term of Φ*M*ð Þ*t* as

$$\begin{split} &\alpha\_{\boldsymbol{\xi}}^{a}P\_{\boldsymbol{\xi},\boldsymbol{\bar{\xi}}}^{(a,b)}(t) = \frac{1}{\Gamma\alpha} \Big|\_{0}^{\boldsymbol{\xi}} (\boldsymbol{\xi}-\boldsymbol{s})^{a-1} P\_{\boldsymbol{\xi},\boldsymbol{\bar{\xi}}}^{(a,b)}(\boldsymbol{\xi}) d\boldsymbol{s} \\ &= \frac{1}{\Gamma\alpha} \Big|\_{0}^{\boldsymbol{\xi}} (\boldsymbol{\xi}-\boldsymbol{s})^{a-1} \sum\_{k=0}^{i} \frac{(-1)^{i-k} (i+b+1) \Gamma(i+k+a+b+1) \boldsymbol{\xi}^{k}}{\Gamma(k+b+1) \Gamma(i+a+b+1) \Gamma(i-k+1) \Gamma(k+1) \boldsymbol{\xi}^{k}} d\boldsymbol{s} \\ &= \frac{1}{\Gamma\alpha} \sum\_{k=0}^{i} \frac{(-1)^{i-k} (i+b+1) \Gamma(i+k+a+b+1)}{\Gamma(k+b+1) \Gamma(i+a+b+1) \Gamma(i-k+1) \Gamma(k+1) \zeta^{k}} \Big|\_{0}^{\boldsymbol{\xi}} (\boldsymbol{\xi}-\boldsymbol{s})^{a-1} \boldsymbol{\xi}^{k} d\boldsymbol{s}. \end{split} \tag{26}$$

Now by applying Laplace transform, we get

$$\begin{split} \mathcal{L}\left(\int\_{0}^{\xi} (\xi - s)^{a-1} s^{k} ds\right) &= \frac{\Gamma(a)\Gamma(k+1)}{s^{k+a+1}} \\ \Rightarrow \int\_{0}^{\xi} (\xi - s)^{a-1} s^{k} ds &= \frac{\Gamma(k+1)\Gamma a}{\Gamma(a+1)} .\end{split} \tag{27}$$

Now putting (27) in (26), we get

$$\,\_0\mathbf{I}\_\xi^q P\_{\xi,i}^{(a,b)}(t) = \sum\_{k=0}^i \frac{(-\mathbf{1})^{i-k} (i+b+1) \Gamma(i+a+b+k+1) \xi^a}{\Gamma(k+b+1) \Gamma(i+a+b+1) \Gamma(i-k+1) \Gamma(k+a)}.\tag{28}$$

Now if *ϕ*ðÞ¼ *t t <sup>n</sup>*, then

$$c\phi(t)\_0\mathbf{I}\_{\xi}^{a}P\_{\xi,i}^{(a,b)}(t) = \sum\_{k=0}^{i} \frac{(-\mathbf{1})^{i-k}(i+b+\mathbf{1})\Gamma(i+a+b+k+\mathbf{1})\xi^{a}ct^{\imath}}{\Gamma(k+b+\mathbf{1})\Gamma(i+a+b+\mathbf{1})\Gamma(i-k+\mathbf{1})\Gamma(k+a)}.\tag{29}$$

Representing *t <sup>n</sup>* in terms of shifted Jacobi polynomials as

$$t^n = \sum\_{j=0}^{M} B\_j P\_{\xi,j}^{(a,b)}(t). \tag{30}$$

By the use of orthogonality relation (11), we have

$$\begin{split} B\_{j} &= \frac{1}{\Omega^{(a,b)}\_{\xi,j}} \Big|\_{0}^{\xi} t^{a} P^{(a,b)}\_{\xi,j}(t) \mathcal{W}^{(a,b)}\_{\xi}(\xi) dt \\ &= \frac{1}{\Omega^{(a,b)}\_{\xi,j}} \sum\_{l=0}^{j} \frac{(-1)^{(\ j-l)} \Gamma(\ j+b+1) \Gamma(\ j+l+a+b+1)}{\Gamma(l+b+1) \Gamma(\ j+a+b+1) \Gamma(\ j-l+1) \Gamma(l+1) \xi^{l}} \Big|\_{0}^{\xi} t^{a+b} (\xi-t)^{a} dt. \end{split} \tag{31}$$

Now by means of Laplace transform, we have

$$\begin{split} \mathcal{L}\left(\int\_{0}^{\xi} t^{a+b} (\xi - t)^{a} dt\right) &= \frac{\Gamma(n+b+1)\Gamma(a+1)}{\xi^{n+b+a+2}} \\ \Rightarrow \quad \int\_{0}^{\xi} t^{n+b} (\xi - t)^{a} dt &= \frac{\Gamma(n+b+1)\Gamma(a+1)}{\Gamma(n+a+b+1)} \xi^{n+a+b+1} .\end{split} \tag{32}$$

Thus putting (32) in (31) and using (12), we get

$$B\_{j} = \sum\_{l=0}^{j} \frac{(-1)^{j-l} (2j+a+b+1) \Gamma(j+1) \Gamma(j+l+a+b+1) \Gamma(n+b+1) \Gamma(a+1) \xi^{n-l}}{\Gamma(j+a+1) \Gamma(l+b+1) \Gamma(j-l+1) \Gamma(l+1) \Gamma(n+a+b+1)}.\tag{33}$$

Now (29) implies that

$$\begin{split} \phi(t)\_0 I\_{\xi}^{a} P\_{\xi,i}^{(a,b)}(t) &= \sum\_{j=0}^{M} \sum\_{k=0}^{i} \frac{(-1)^{i-k} (i+b+1) \Gamma(i+a+b+k+1) \xi^a c}{\Gamma(k+b+1) \Gamma(i+a+b+1) \Gamma(i-k+1) \Gamma(k+a)} B\_{j} P\_{\xi,j}^{(a,b)}(t) \\ &= \sum\_{j=0}^{M} \Delta\_{i,k,a} B\_{j} P\_{\xi,j}^{(a,b)}(t) = \sum\_{j=0}^{M} \Upsilon\_{i,j,k} P\_{\xi,j}^{(a,b)}(t), \end{split} \tag{34}$$

where

$$\Upsilon\_{i,j,k} = \sum\_{j=0}^{M} \Delta\_{i,k,a} B\_j,$$

$$\Delta\_{i,k,a} = \sum\_{k=0}^{i} \frac{(-1)^{i-k} (i+b+1) \Gamma(i+a+b+k+1) \xi^a c}{\Gamma(k+b+1) \Gamma(i+a+b+1) \Gamma(i+k+1) \Gamma(k+a)}.$$

Evaluating the result (34) for different values of *j* we get the required result.

#### **4. Applications of operational matrices**

In this section, we show fundamental importance of operational matrices of fractional order derivative and integration. We apply them to solve some multiterms fractional order boundary value problems of fractional differential equations. Consider the following general fractional differential equation with constant coefficient and given boundary conditions

$$\mathbf{D}^{\alpha}\mathbf{x}(t) + A\mathbf{D}^{\nu}\mathbf{x}(t) + B\mathbf{D}^{\alpha}\mathbf{x}(t) + \mathbf{C}\mathbf{x}(t) = f(t), t \in [0, T], A, B, \mathbf{C} \in \mathbb{R},\tag{35}$$

subject the boundary conditions

$$\mathfrak{x}(\mathbf{0}) = \mathfrak{x}\_0, \quad \mathfrak{x}(T) = \mathfrak{x}\_1,$$

where *f t*ð Þ∈*C*½ � 0, *T* is a source terms and 1<*α* ≤2, 0<*ν*,*ω* <1*:*.

To obtain the solutions of the (35) in terms of shifted Jacobi polynomials, we proceed as

$$D^a \mathfrak{x}(t) = \mathbf{K}\_M^T \Phi\_M(t). \tag{36}$$

*Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order DOI: http://dx.doi.org/10.5772/intechopen.102054*

Applying I*<sup>α</sup>*, definition (7), and lemma (2.2.1), then we get

$$\mathbf{x}(t) = \mathbf{K}\_{\mathcal{M}}^T \mathbf{H}\_{\mathcal{M} \times \mathcal{M}}^{(a)} \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \boldsymbol{c}\_0 + \boldsymbol{c}\_1 t. \tag{37}$$

By means of boundary conditions we get

$$\mathbf{x}(t) = \mathbf{K}\_{M}^{T} \mathbf{H}\_{M \times M}^{(a)} \boldsymbol{\Phi}\_{M}(t) + \boldsymbol{x}\_{0} + t \frac{(\boldsymbol{\varkappa}\_{1} - \boldsymbol{\varkappa}\_{0})}{T} - \frac{t}{T} \mathbf{K}\_{M}^{T} \mathbf{H}\_{M \times M}^{(a)} \boldsymbol{\Phi}\_{M}(t) \Big|\_{t=T} \tag{38}$$

using the approximation *x*<sup>0</sup> þ *t* ð Þ *x*1�*x*<sup>0</sup> *<sup>T</sup>* <sup>≃</sup>F<sup>1</sup> *<sup>M</sup>*Φ*M*ð Þ*<sup>t</sup>* , *<sup>t</sup> <sup>T</sup>* ¼ *ϕ*ð Þ*t* , also using *<sup>c</sup>ϕ*ð Þ*<sup>t</sup>* <sup>K</sup>*<sup>T</sup> M*Hð Þ *<sup>α</sup> <sup>M</sup>*�*M*Φ*M*ð Þ*<sup>t</sup>* � � � *<sup>t</sup>*¼*<sup>T</sup>* <sup>¼</sup> <sup>K</sup>*<sup>T</sup> M*Qð Þ *<sup>α</sup>*,*<sup>ϕ</sup> <sup>M</sup>*�*M*Φ*M*ð Þ*<sup>t</sup>* , then (38) implies that

$$\begin{split} \mathbf{x}(t) &= \mathbf{K}\_{\mathcal{M}}^{T} \mathbf{H}\_{\mathcal{M}\times\mathcal{M}}^{(a)} \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{\mathcal{M}}^{1} \boldsymbol{\Phi}(t) - c\boldsymbol{\phi}(t) \mathbf{K}\_{\mathcal{M}}^{T} \mathbf{H}\_{\mathcal{M}\times\mathcal{M}}^{(a)} \boldsymbol{\Phi}\_{\mathcal{M}}(t) \Big|\_{t=T} \\ &= \mathbf{K}\_{\mathcal{M}}^{T} \mathbf{H}\_{\mathcal{M}\times\mathcal{M}}^{(a)} \boldsymbol{\Phi}\_{\mathcal{M}}(t) - \mathbf{K}\_{\mathcal{M}}^{T} \mathbf{Q}\_{\mathcal{M}\times\mathcal{M}}^{(a,\phi)} \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{\mathcal{M}}^{1} \boldsymbol{\Phi}\_{\mathcal{M}}(t) \\ &= \mathbf{K}\_{\mathcal{M}}^{T} \left( \mathbf{H}\_{\mathcal{M}\times\mathcal{M}}^{a} - \mathbf{Q}\_{\mathcal{M}\times\mathcal{M}}^{(a,\phi)} \right) \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{\mathcal{M}}^{1} \boldsymbol{\Phi}\_{\mathcal{M}}(t). \end{split} \tag{39}$$

Now from (39) we have

$$\mathbf{D}^{\nu}\mathbf{x}(t) = \mathbf{K}\_{M}^{T} \left( \mathbf{H}\_{M\times M}^{(a)} - \mathbf{Q}\_{M\times M}^{(a,\phi)} \right) \mathbf{G}\_{M\times M}^{(\nu)} \boldsymbol{\Phi}\_{M}(t) + \mathbf{F}\_{M}^{1} \mathbf{G}\_{M\times M}^{(\nu)} \boldsymbol{\Phi}\_{M}(t)$$
 
$$\mathbf{D}^{\alpha}\mathbf{x}(t) = \mathbf{K}\_{M}^{T} \left( \mathbf{H}\_{M\times M}^{(a)} - \mathbf{Q}\_{M\times M}^{(a,\phi)} \right) \mathbf{G}\_{M\times M}^{(\alpha)} \boldsymbol{\Phi}\_{M}(t) + \mathbf{F}\_{M}^{1} \mathbf{G}\_{M\times M}^{(\alpha)} \boldsymbol{\Phi}\_{M}(t) \tag{40}$$
 
$$\text{and approximating the source term as} \qquad f(t) \simeq \mathbf{F}\_{M}^{2} \boldsymbol{\Phi}\_{M}(t).$$

Putting (36)–(40) in (35), which yields

K*T <sup>M</sup>*Φ*M*ð Þþ*<sup>t</sup> <sup>A</sup>* <sup>K</sup>*<sup>T</sup> <sup>M</sup>* <sup>H</sup>ð Þ *<sup>α</sup> <sup>M</sup>*�*<sup>M</sup>* � <sup>Q</sup>ð Þ *<sup>α</sup>*,*<sup>ϕ</sup> M*�*M* � �Gð Þ*<sup>ν</sup> <sup>M</sup>*�*<sup>M</sup>*Φ*M*ðÞþ*<sup>t</sup>* <sup>F</sup><sup>1</sup> *M*Gð Þ*<sup>ν</sup> <sup>M</sup>*�*<sup>M</sup>*Φ*M*ð Þ*<sup>t</sup>* h i <sup>þ</sup>*<sup>B</sup>* <sup>K</sup>*<sup>T</sup> <sup>M</sup>* <sup>H</sup>ð Þ *<sup>α</sup> <sup>M</sup>*�*<sup>M</sup>* � <sup>Q</sup>ð Þ *<sup>α</sup>*,*<sup>ϕ</sup> M*�*M* � �Gð Þ *<sup>ω</sup> <sup>M</sup>*�*<sup>M</sup>*Φ*M*ðÞþ*<sup>t</sup>* <sup>F</sup><sup>1</sup> *M*Gð Þ *<sup>ω</sup> <sup>M</sup>*�*<sup>M</sup>*Φ*M*ð Þ*<sup>t</sup>* h i <sup>þ</sup>*<sup>C</sup>* <sup>K</sup>*<sup>T</sup> <sup>M</sup>* <sup>H</sup>ð Þ *<sup>α</sup> <sup>M</sup>*�*<sup>M</sup>* � <sup>Q</sup>ð Þ *<sup>α</sup>*,*<sup>ϕ</sup> M*�*M* � �Φ*M*ðÞþ*<sup>t</sup>* <sup>F</sup><sup>1</sup> *<sup>M</sup>*Φ*M*ð Þ*t* h i � F2 *<sup>M</sup>*Φ*M*ðÞ¼ *t* 0 K*T <sup>M</sup>*Φ*M*ð Þþ*<sup>t</sup>* <sup>K</sup>*<sup>T</sup> <sup>M</sup>* <sup>H</sup>ð Þ *<sup>α</sup> <sup>M</sup>*�*<sup>M</sup>* � <sup>Q</sup>ð Þ *<sup>α</sup>*,*<sup>ϕ</sup> M*�*M* � � *<sup>A</sup>*Gð Þ*<sup>ν</sup> <sup>M</sup>*�*<sup>M</sup>* <sup>þ</sup> *<sup>B</sup>*Gð Þ *<sup>ω</sup> <sup>M</sup>*�*<sup>M</sup>* <sup>þ</sup> *CI* h iΦ*M*ð Þ*<sup>t</sup>* <sup>þ</sup> F1 *<sup>M</sup> <sup>A</sup>*Gð Þ*<sup>ν</sup> <sup>M</sup>*�*<sup>M</sup>* <sup>þ</sup> *<sup>B</sup>*Gð Þ *<sup>ω</sup> M*�*M* � � <sup>þ</sup> *<sup>C</sup>*F1 *<sup>M</sup>* � F2 *M* h iΦ*M*ðÞ¼ *<sup>t</sup>* <sup>0</sup> K*T <sup>M</sup>* <sup>þ</sup> <sup>K</sup>*<sup>T</sup> <sup>M</sup>* H*<sup>α</sup> <sup>M</sup>*�*<sup>M</sup>* � <sup>Q</sup>ð Þ *<sup>α</sup>*,*<sup>ϕ</sup> M*�*M* � � *<sup>A</sup>*Gð Þ*<sup>ν</sup> <sup>M</sup>*�*<sup>M</sup>* <sup>þ</sup> *<sup>B</sup>*Gð Þ *<sup>ω</sup> <sup>M</sup>*�*<sup>M</sup>* <sup>þ</sup> *CI* h i <sup>þ</sup>F<sup>1</sup> *<sup>M</sup> <sup>A</sup>*Gð Þ*<sup>ν</sup> <sup>M</sup>*�*<sup>M</sup>* <sup>þ</sup> *<sup>B</sup>*Gð Þ *<sup>ω</sup> <sup>M</sup>*�*<sup>M</sup>* <sup>þ</sup> *CI* � � � F2 *<sup>M</sup>* ¼ 0*:* (41)

Which is a simple algebraic equation. Solving this equation for **K***<sup>M</sup>* and putting it in (39) we get the required approximate solution of the corresponding boundary value problem.

#### **4.1 Coupled system of boundary value problems for fractional order differential equations**

In this subsection, we use operational matrices to derive procedure for the numerical solutions of coupled system. We consider the following general coupled system of FDEs as

$$\begin{aligned} \mathbf{D}^{\alpha}\mathbf{x}(t) + A\_1 \mathbf{D}^{\nu\_1}\mathbf{x}(t) + B\_1 \mathbf{D}^{\nu\_2}\mathbf{y}(t) + C\_1 \mathbf{x}(t) + D\_1 \mathbf{y}(t) &= f(t), \ t \in [0, T], \\ \mathbf{D}^{\beta}\mathbf{y}(t) + A\_2 \mathbf{D}^{\alpha\_1}\mathbf{x}(t) + B\_2 \mathbf{D}^{\alpha\_2}\mathbf{y}(t) + C\_2 \mathbf{x}(t) + D\_2 \mathbf{y}(t) &= \mathbf{g}(t), \ t \in [0, T], \end{aligned} \tag{42}$$

subject to the boundary conditions given by

$$\begin{aligned} \varkappa(\mathbf{0}) &= \varkappa\_0, \ \varkappa(T) = \varkappa\_1, \\ \varkappa(\mathbf{0}) &= \mathcal{y}\_0, \ \mathcal{y}(T) = \mathcal{y}\_1, \end{aligned} \tag{43}$$

where 1<*α*, *β* ≤2, 0< *νi*,*ω<sup>i</sup>* <1ð Þ *i* ¼ 1, 2 *:* Further, *f t*ð Þ, *g t*ð Þ∈*C*½ � 0, *T* are the source term of the system (42), (43) and *x*0, *y*0, *x*1, *y*<sup>1</sup> are real constants, *Ai*, *Bi*,*Ci*, *Di*ð Þ *i* ¼ 1, 2 are any constants. To approximate the above system in terms of shifted Jacobi polynomials, we take

$$\mathbf{D}^a \mathbf{x}(t) = \mathbf{K}\_M^T \boldsymbol{\Phi}\_M(t) \tag{44}$$

and

$$\mathbf{D}^{\beta}y(t) = \mathbf{L}\_{\mathcal{M}}^{T}\boldsymbol{\Phi}\_{\mathcal{M}}(t) \tag{45}$$

Applying I*<sup>α</sup>*, I*<sup>β</sup>* and lemma (1) and boundary conditions, the above Eqs. (44), (45) imply that

$$\mathbf{x}(t) = \mathbf{K}\_{M}^{T} \left( \mathbf{H}\_{M \times M}^{(a)} - \mathbf{Q}\_{M \times M}^{(a, \phi)} \right) \Phi\_{M}(t) + \mathbf{F}\_{M}^{1} \Phi\_{M}(t) \tag{46}$$

$$\mathbf{y}(t) = \mathbf{L}T\_M \Big(\mathbf{H}\_{M \times M}^{(\beta)} - \mathbf{Q}\_{M \times M}^{(\beta, \phi)}\Big) \Phi\_M(t) + \mathbf{F}\_M^2 \Phi\_M(t). \tag{47}$$

Now taking fractional order derivative of (46), (47), we have

$$\begin{aligned} \mathbf{D}^{\boldsymbol{\nu}\_{1}}\boldsymbol{\mathbf{x}}(t) &= \mathbf{K}\_{\mathcal{M}}^{T} \Big( \mathbf{H}\_{\mathcal{M}\times\mathcal{M}}^{(a)} - \mathbf{Q}\_{\mathcal{M}\times\mathcal{M}}^{(a,\phi)} \Big) \mathbf{G}\_{\mathcal{M}\times\mathcal{M}}^{(\nu\_{1})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{\mathcal{M}}^{1} \mathbf{G}\_{\mathcal{M}\times\mathcal{M}}^{(\nu\_{1})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) \\ \mathbf{D}^{\boldsymbol{\alpha}\boldsymbol{\eta}}\boldsymbol{x}(t) &= \mathbf{K}\_{\mathcal{M}}^{T} \Big( \mathbf{H}\_{\mathcal{M}\times\mathcal{M}}^{(a)} - \mathbf{Q}\_{\mathcal{M}\times\mathcal{M}}^{(a,\phi)} \Big) \mathbf{G}\_{\mathcal{M}\times\mathcal{M}}^{(\boldsymbol{\alpha}\boldsymbol{\eta})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{\mathcal{M}}^{1} \mathbf{G}\_{\mathcal{M}\times\mathcal{M}}^{(\boldsymbol{\alpha}\boldsymbol{\eta})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) \end{aligned} \tag{48}$$

and

$$\begin{split} \mathbf{D}^{\nu\_{2}}\boldsymbol{\mathcal{y}}(t) &= \mathbf{L}\_{\mathcal{M}}^{T} \left( \mathbf{H}\_{\mathcal{M}\times\mathcal{M}}^{(\boldsymbol{\beta})} - \mathbf{Q}\_{\mathcal{M}\times\mathcal{M}}^{(\boldsymbol{\beta},\boldsymbol{\phi})} \right) \mathbf{G}\_{\mathcal{M}\times\mathcal{M}}^{(\nu\_{2})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{\mathcal{M}}^{2} \mathbf{G}\_{\mathcal{M}\times\mathcal{M}}^{(\nu\_{2})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) \\ \mathbf{D}^{(\boldsymbol{\alpha}\_{2})}\boldsymbol{\jmath}(t) &= \mathbf{L}\_{\mathcal{M}}^{T} \left( \mathbf{H}\_{\mathcal{M}\times\mathcal{M}}^{(\boldsymbol{\beta})} - \mathbf{Q}\_{\mathcal{M}\times\mathcal{M}}^{(\boldsymbol{\beta},\boldsymbol{\phi})} \right) \mathbf{G}\_{\mathcal{M}\times\mathcal{M}}^{(\boldsymbol{\alpha}\_{2})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{\mathcal{M}}^{2} \mathbf{G}\_{\mathcal{M}\times\mathcal{M}}^{(\boldsymbol{\alpha}\_{2})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) . \end{split} \tag{49}$$

Putting (45)–(49) together with the use of approximation *f t*ð Þ<sup>≃</sup> F3 *<sup>M</sup>*Φ*M*ð Þ*<sup>t</sup>* , *g t*ð Þ≃F<sup>4</sup> *<sup>M</sup>*Φ*M*ð Þ*t* in (42), we obtain

$$\begin{aligned} \mathbf{K}\_{M}^{T}\boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{A}\_{1} \Big[ \mathbf{K}\_{M}^{T} \Big( \mathbf{H}\_{M \times M}^{(a)} - \mathbf{Q}\_{M \times M}^{(a,\phi)} \Big) \mathbf{G}\_{M \times M}^{(\nu\_{1})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{M}^{1} \mathbf{G}\_{M \times M}^{(\nu\_{1})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) \Big] \\ + \mathbf{B}\_{1} \Big[ \mathbf{L}\_{M}^{T} \Big( \mathbf{H}\_{M \times M}^{(\beta)} - \mathbf{Q}\_{M \times M}^{(\beta,\phi)} \Big) \mathbf{G}\_{M \times M}^{(\nu\_{2})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{M}^{2} \mathbf{G}\_{M \times M}^{(\nu\_{2})} \boldsymbol{\Phi}\_{\mathcal{M}}(t) \Big] \\ + \mathbf{C}\_{1} \Big[ \mathbf{K}\_{M}^{T} \Big( \mathbf{H}\_{M \times M}^{(a)} - \mathbf{Q}\_{M \times M}^{(a,\phi)} \Big) \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{M}^{1} \boldsymbol{\Phi}\_{\mathcal{M}}(t) \Big] \\ + \mathbf{D}\_{1} \Big[ \mathbf{L}\_{M}^{T} \Big( \mathbf{H}\_{M \times M}^{(\beta)} - \mathbf{Q}\_{M \times M}^{(\beta,\phi)} \Big) \boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{F}\_{M}^{2} \boldsymbol{\Phi}\_{\mathcal{M}}(t) \Big] \\ - \mathbf{F}\_{M}^{3} \boldsymbol{\Phi}\_{\mathcal{M}}(t) = \mathbf{0} \end{aligned} \tag{50}$$

*Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order DOI: http://dx.doi.org/10.5772/intechopen.102054*

$$\begin{split} \mathbf{L}\_{M}^{T}\boldsymbol{\Phi}\_{\mathcal{M}}(t) + \mathbf{A}\_{2} \left[ \mathbf{K}\_{M}^{T} \left( \mathbf{H}\_{M\times M}^{(a)} - \mathbf{Q}\_{M\times M}^{(a,\phi)} \right) \mathbf{G}\_{M\times M}^{(a\alpha)} \boldsymbol{\Phi}\_{M}(t) + \mathbf{F}\_{M}^{1} \mathbf{G}\_{M\times M}^{(a\alpha)} \boldsymbol{\Phi}\_{M}(t) \right] \\ + \mathbf{B}\_{2} \left[ \mathbf{L}\_{M}^{T} \left( \mathbf{H}\_{M\times M}^{(\beta)} - \mathbf{Q}\_{M\times M}^{(\beta,\phi)} \right) \mathbf{G}\_{M\times M}^{(\alpha\_{1})} \boldsymbol{\Phi}\_{M}(t) + \mathbf{F}\_{M}^{2} \mathbf{G}\_{M\times M}^{(\alpha\_{2})} \boldsymbol{\Phi}\_{M}(t) \right] \\ + \mathbf{C}\_{2} \left[ \mathbf{K}\_{M}^{T} \left( \mathbf{H}\_{M\times M}^{(a)} - \mathbf{Q}\_{M\times M}^{(a,\phi)} \right) \boldsymbol{\Phi}\_{M}(t) + \mathbf{F}\_{M}^{1} \boldsymbol{\Phi}\_{M}(t) \right] \\ + \mathbf{D}\_{2} \left[ \mathbf{L}\_{M}^{T} \left( \mathbf{H}\_{M\times M}^{(\beta)} - \mathbf{Q}\_{M\times M}^{(\beta,\phi)} \right) \boldsymbol{\Phi}\_{M}(t) + \mathbf{F}\_{M}^{2} \boldsymbol{\Phi}\_{M}(t) \right] \\ - \mathbf{F}\_{M}^{4} \boldsymbol{\Phi}\_{M}(t) = 0. \end{split}$$

Rearranging the terms in above system (50) and using the notation for simplicity

$$\begin{aligned} \mathbf{S} &= \mathbf{F}\_M^1 \left( A\_1 \mathbf{G}^{(\nu\_1)} + \mathbf{C}\_1 I \right) + \mathbf{F}\_M^2 \left( B\_1 \mathbf{G}^{(\nu\_2)} + D\_1 I \right) - \mathbf{F}\_M^3 \\ \mathbf{R} &= \mathbf{F}\_M^1 \left( A\_2 \mathbf{G}^{(\nu\_1)} + \mathbf{C}\_2 I \right) + \mathbf{F}\_M^2 \left( B\_2 \mathbf{G}^{(\nu\_2)} + D\_2 I \right) - \mathbf{F}\_M^4 \\ \mathbf{H}^{(a)} - \mathbf{Q}^{(a, \phi)} &= \widehat{\mathbf{P}} \quad , \quad \mathbf{H}^{(\beta)} - \mathbf{Q}^{(\beta, \phi)} = \widehat{\mathbf{Q}} \end{aligned} \tag{50a}$$

then witting in matrix form, we have

K*T <sup>M</sup>* L*<sup>T</sup> M* � � Φ*M*ð Þ*t* 0 0 Φ*M*ð Þ*t* " # <sup>þ</sup> <sup>K</sup>*<sup>T</sup> <sup>M</sup>* L*<sup>T</sup> M* � � <sup>P</sup> z}|{ *<sup>A</sup>*1G ð Þ *<sup>ν</sup>*<sup>1</sup> <sup>þ</sup> *<sup>C</sup>*1*<sup>I</sup>* <sup>Q</sup> z}|{ *<sup>B</sup>*1G ð Þ *<sup>ν</sup>*<sup>2</sup> <sup>þ</sup> *<sup>D</sup>*1*<sup>I</sup>* P z}|{ *<sup>A</sup>*2G ð Þ *<sup>ω</sup>*<sup>1</sup> <sup>þ</sup> *<sup>C</sup>*2*<sup>I</sup>* <sup>Q</sup> z}|{ *<sup>B</sup>*2G ð Þ *<sup>ω</sup>*<sup>2</sup> <sup>þ</sup> *<sup>C</sup>*2*<sup>I</sup>* 2 6 4 3 7 5 Φ*M*ð Þ*t* 0 0 Φ*M*ð Þ*t* " # þ½ � S R Φ*M*ð Þ*t* 0 0 Φ*M*ð Þ*t* " # ¼ 0 ) <sup>K</sup>*<sup>T</sup> <sup>M</sup>* L*<sup>T</sup> M* � � <sup>þ</sup> <sup>K</sup>*<sup>T</sup> <sup>M</sup>* L*<sup>T</sup> M* � � <sup>P</sup> z}|{ *<sup>A</sup>*1G ð Þ *<sup>ν</sup>*<sup>1</sup> <sup>þ</sup> *<sup>C</sup>*1*<sup>I</sup>* <sup>Q</sup> z}|{ *<sup>B</sup>*1G ð Þ *<sup>ν</sup>*<sup>2</sup> <sup>þ</sup> *<sup>D</sup>*1*<sup>I</sup>* P z}|{ *<sup>A</sup>*2G ð Þ *<sup>ω</sup>*<sup>1</sup> <sup>þ</sup> *<sup>C</sup>*2*<sup>I</sup>* <sup>Q</sup> z}|{ *<sup>B</sup>*2G ð Þ *<sup>ω</sup>*<sup>2</sup> <sup>þ</sup> *<sup>C</sup>*2*<sup>I</sup>* 2 6 4 3 7 <sup>5</sup> <sup>þ</sup> ½ �¼ S R <sup>0</sup> (51)

which is an algebraic equation and can easily be solved for unknown matrices K*M*, L*M*, then putting it in (38) and (39), we get required approximations for *x t*ð Þ, *y t*ð Þ.

#### **5. Numerical examples**

In this section, we will apply our proposed scheme to some practical problems. **Example 1.** *Consider the problems of Buckling of a thin vertical column such that when a constant compressive force F is applied to the said thin column of length L hanged at both ends whose general differential equations in the absence of source term are given by*

$$\begin{aligned} \text{EL}\frac{d^a y}{d\mathbf{x}^a} + \text{Fy} &= f(t), \ \mathbf{1} < a \le 2, \\ \text{subject to the boundary conditions} \quad &y(\mathbf{0}) = y\_1, \ y(L) = y\_2, \end{aligned} \tag{52}$$


**Table 1.**

*Maximum absolute error for T* ¼ 1, *α* ¼ 1*:*8, *ν*<sup>1</sup> ¼ *ν*<sup>2</sup> ¼ 0 *for example, (1).*

where *E*,*I* are Young modulus and moment of inertia, respectively. Where *E* ¼ <sup>2</sup>*:*<sup>6</sup> � <sup>10</sup><sup>7</sup>*lb=in*,*<sup>I</sup>* <sup>¼</sup> <sup>0</sup>*:*25*π*<sup>2</sup>*R*4, where *<sup>R</sup>* is the radius of the column. Taking *<sup>L</sup>* <sup>¼</sup> 1*in*, *R* ¼ 1*in* the exact deflection at *α* ¼ 2 is *y x*ð Þ¼ sin ð Þ *πx .* We apply our technique introduced in section 4ð Þ for various choices of *a*, *b*, *M.* The observations are given in **Table 1** for the absolute error ∣*yapp* � *yexact*∣*.*

In **Table 2**, we have shown the maximum absolute error for various choices of *a*, *b*, *M.* Also in **Figures 1** and **2**, we have shown the images of their comparison between exact and approximate values. Further, we observed that as the *α* ! 2, the error is reducing and one can check that at *α* ¼ 2, the absolute error is below 1*:*8 � <sup>10</sup>�<sup>13</sup> as given in the next table. From above tables, we observe that as order *<sup>α</sup>* ! <sup>2</sup>


#### **Table 2.** *Maximum absolute error at a* ¼ *b* ¼ 0, *and for different α*, *ν*1, *ν*<sup>2</sup> *for example, (2).*

#### **Figure 1.**

*Subplot (a) represents comparison of approximate and exact y x*ð Þ *of example 1. Stars and dots represent exact solution while the red curve represents approximate solution. Setting α* ¼ 1*:*8, *M* ¼ 4, *T* ¼ 1, *a* ¼ *b* ¼ 1*: the subplot represents their absolute error at M* ¼ 8.

*Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order DOI: http://dx.doi.org/10.5772/intechopen.102054*

**Figure 2.**

*Subplot (a) represents comparison of approximate and exact y x*ð Þ *of example 1.Red dots represent exact solution while the blue curve represents approximate solution. Setting α* ¼ 1*:*8, *a* ¼ *b* ¼ 0, *M* ¼ 8, *T* ¼ 1*: the subplot (b) represents their absolute error at M* ¼ 8.

the approximation converges to the value at *α* ¼ 2*.* This phenomenon shows the accuracy of the approximate solutions.

**Example 2.** *Consider the following general problem*

$$\begin{aligned} \text{D}^{a}Y + A\text{D}^{\nu\_{1}}Y + B\text{D}^{\nu\_{2}}Y + CY &= f(\mathbf{x}), \ \mathbf{1} < a \le \mathbf{2}, \ \mathbf{0} < \nu\_{1}, \nu\_{2} \le \mathbf{1}, \\\text{subject to the boundary conditions} \qquad \mathbf{Y}(\mathbf{0}) = \mathbf{1}, \ \mathbf{Y}(\mathbf{1}) = \mathbf{5}. \end{aligned} \tag{53}$$

Let the exact solution *Y x*ð Þ¼ *<sup>x</sup>*<sup>4</sup> <sup>þ</sup> *<sup>x</sup>*<sup>3</sup> <sup>þ</sup> *<sup>x</sup>*<sup>2</sup> <sup>þ</sup> *<sup>x</sup>* <sup>þ</sup> 1, and *<sup>A</sup>* <sup>¼</sup> 1, *<sup>B</sup>* <sup>¼</sup> 2,*<sup>C</sup>* <sup>¼</sup> 3, the source terms

$$\begin{split} f(\mathbf{x}) &= 3\mathbf{x}^4 + 3\mathbf{x}^3 + 15\mathbf{x}^2 + 9\mathbf{x} + 5 + \frac{24}{\Gamma} \frac{\mathbf{x}^4}{5} + \frac{6}{\Gamma} \frac{\mathbf{x}^4}{5} + \frac{2}{\Gamma} \frac{\mathbf{x}^4}{5} + \frac{1}{\Gamma} \frac{\mathbf{x}^4}{5} \\ &+ \frac{48}{\Gamma} \frac{23}{5} \mathbf{x}^{\frac{3}{2}} + \frac{12}{\Gamma} \frac{\mathbf{x}^{\frac{3}{2}}}{5} + \frac{4}{\Gamma} \frac{\mathbf{x}^{\frac{3}{2}}}{5} \mathbf{x}^{\frac{3}{2}} .\end{split} \tag{54}$$

We approximate the solutions of this problem with proposed method by setting *α* ¼ 2, *ν*<sup>1</sup> ¼ 0*:*8, *ν*<sup>2</sup> ¼ 0*:*4, *a* ¼ *b* ¼ 1, *M* ¼ 8, *T* ¼ 1*:* and observed that the scheme provides much more accurate approximations for very small-scale level even at scale level *<sup>M</sup>* <sup>¼</sup> 8 maximum absolute error is 9 � <sup>10</sup>�<sup>17</sup> as shown in **Figure 3**.

The absolute error for different choices of *a*, *b* and *M* and for *α* ¼ 2, *ν*<sup>1</sup> ¼ 0*:*8, *ν*<sup>2</sup> ¼ 0*:*4 is given in **Table 3**.

From **Table 4**, we observed that as the orders of derivatives tend to their corresponding integer values, the approximate solutions tend to its exact value, which demonstrate the accuracy of numerical solutions obtained in our proposed method.

In the following examples, we solve some coupled systems under the given conditions by our proposed method.

**Example 3.** Consider the coupled system given by

$$\begin{aligned} \mathbf{D}^{\alpha}X &+ A\_1 \mathbf{D}^{\alpha}X + B\_1 \mathbf{D}^{\alpha}Y + C\_1 X + D\_1 Y = f(\mathbf{x}), \ 1 & 0 < \alpha \le 2, \ 0 < \nu\_1, \nu\_2 \le 1, \\\mathbf{D}^{\beta}Y &+ A\_2 \mathbf{D}^{\alpha}X + B\_2 \mathbf{D}^{\alpha}Y + C\_2 X + D\_2 Y = \mathbf{g}(\mathbf{x}), \ 1 & 0 \le \mathbf{2}, \ 0 < \nu\_1, \nu\_2 \le 1, \end{aligned} \tag{55}$$

#### **Figure 3.**

*Subplot (a) represents comparison of approximate and exact Y x*ð Þ *of example 2. Blue curve represents exact solution while the red square dots represent approximate solution. The subplot (b) represents their absolute error at M* ¼ 8.


**Table 3.**

*Maximum absolute error at various values of α*, *a*, *b and M for example, (1).*


**Table 4.** *Maximum absolute error at α* ¼ 2, *ν*<sup>1</sup> ¼ 0*:*8, *ν*<sup>2</sup> ¼ 0*:*4 *for example, (2).*

subject to the boundary conditions

$$\begin{aligned} X(\mathbf{0}) &= \mathbf{0}.\mathbf{5}, & X(\mathbf{1}) &= \mathbf{0}.\mathbf{5} \\ Y(\mathbf{0}) &= \mathbf{0}.\mathbf{6}, & Y(\mathbf{1}) &= \mathbf{0}.\mathbf{6}.\end{aligned} \tag{56}$$

We solve this problem under the given parameters

$$\begin{aligned} A\_1 &= 2, A\_2 = 2, B\_1 = 2, B\_2 = 2, C\_1 = 3, C\_2 = 4, D\_1 = 5, D\_2 = 5, \\\ a &= \beta = 2, \nu\_1 = a\_2 = 1, \nu\_1 = a\_1 = 1, a = b = 0. \end{aligned}$$

*Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order DOI: http://dx.doi.org/10.5772/intechopen.102054*

Let the exact solutions at *α* ¼ *β* ¼ 2, *ν<sup>i</sup>* ¼ *ω<sup>i</sup>* ¼ 1ð Þ *i* ¼ 1, 2 are given by

$$X(\mathbf{x}) = \mathbf{x}^2 (T - \mathbf{x})^2 + \mathbf{0.5},\\ Y(\mathbf{x}) = \mathbf{x}(T - \mathbf{x})^3 + \mathbf{0.6},\tag{57}$$

*the source terms are given by*

$$\begin{aligned} f(\mathbf{x}) &= 4\mathbf{x}(2\mathbf{x}-2) - 6\mathbf{x}(\mathbf{x}-\mathbf{1})^2 - 5\mathbf{x}(\mathbf{x}-\mathbf{1})^3 + 2(\mathbf{x}-\mathbf{1})^2 - 2(\mathbf{x}-\mathbf{1})^3 + 3\mathbf{x}^2(\mathbf{x}-\mathbf{1})^2 \\ &+ \frac{\left(3172993274221175\mathbf{x}^{(\frac{5}{3})}(50\mathbf{x}^2 - 85\mathbf{x} + 34)\right)}{3349552228568064} + 2\mathbf{x}^2 + \frac{9}{2} \\ g(\mathbf{x}) &= 2\mathbf{x}^2(2\mathbf{x}-2) - 2\mathbf{x}(\mathbf{x}-\mathbf{1})^2 - 5\mathbf{x}(\mathbf{x}-\mathbf{1})^3 - 6(\mathbf{x}-\mathbf{1})^2 - 2(\mathbf{x}-\mathbf{1})^3 - 3\mathbf{x}(2\mathbf{x}-2) \\ &+ 3\mathbf{x}^2(\mathbf{x}-\mathbf{1})^2 + \frac{9}{2} \end{aligned} \tag{59}$$

We approximate the coupled systems at scale level *M* ¼ 6 for the given parameters. We see from **Figures 4** and **5** that our method works very well, and the

#### **Figure 4.**

*Subplot (a) represents comparison of approximate and exact X x*ð Þ *of example 3. Blue curve represents exact solution while the red circles represent approximate solution. The subplot (b) represents their absolute error at M* ¼ 6.

#### **Figure 5.**

*Subplot (c) represents comparison of approximate and exact Y x*ð Þ *of example 3. Blue dots represent exact solution while the red curve represents approximate solution. The subplot (d) represents their absolute error at M* ¼ 6.


<sup>6</sup> <sup>3</sup> � <sup>10</sup>�<sup>16</sup> <sup>7</sup> � <sup>10</sup>�<sup>16</sup> 6 11 <sup>4</sup>*:*<sup>5</sup> � <sup>10</sup>�<sup>16</sup> <sup>7</sup> � <sup>10</sup>�<sup>16</sup> 8 0*:*5 1 <sup>4</sup>*:*<sup>5</sup> � <sup>10</sup>�<sup>16</sup> <sup>6</sup> � <sup>10</sup>�<sup>16</sup> 8 10*:*<sup>5</sup> <sup>3</sup>*:*<sup>5</sup> � <sup>10</sup>�<sup>16</sup> <sup>1</sup>*:*<sup>2</sup> � <sup>10</sup>�<sup>15</sup>

**Table 5.**

*Recent Advances in Polynomials*

*Maximum absolute error at various values of a*, *b and M for example, (3).*


**Table 6.**

*Maximum absolute error at specific values of a* ¼ *b* ¼ 0 *for example, (3).*

absolute error is about 3 � <sup>10</sup>�<sup>16</sup> for *X x*ð Þ and for *Y x*ð Þ this amount is 7 � <sup>10</sup>�<sup>16</sup>*:* Which is very small numbers, which prove the applicability of our methods. Absolute error for various *a*, *b*, *M* is given in next **Table 5**.

Absolute error for specific *a*, *b* and for various *α*, *β*, *ν*1, *ν*2,*ω*1,*ω*<sup>2</sup> is given in **Table 6**.

From **Table 5**, we see the effect of *a*, *b* at various scale level *M*. While **Table 6** indicates the effect of various choices of *α*, *β*, *ν*1, *ν*2,*ω*1,*ω*<sup>2</sup> at different scale level. As the scale level increases, values of *α*, *β*, *ν*1, *ν*2,*ω*1,*ω*<sup>2</sup> tend to their integer values. The absolute error is reducing and the solutions are tending to their exact value, which demonstrate the applicability of the proposed method.

**Example 4.** Consider the mathematical models of fractionally damped coupled system of spring masses system whose model is given by

$$\begin{aligned} \mathbf{D}^{\theta}X + \frac{\lambda}{m\_{1}}\mathbf{D}^{\nu\_{1}}X - \frac{\mu}{m\_{1}}\mathbf{D}^{\nu\_{2}}Y + \frac{\kappa\_{1}}{m\_{1}}X - \frac{\kappa\_{2}}{m\_{1}}Y &= f(\mathbf{x}), \ \mathbf{1} < a \le 2, \ \mathbf{0} < \nu\_{1}, \nu\_{2} \le 1, \\\ \mathbf{D}^{\beta}Y - \frac{\lambda}{m\_{2}}\mathbf{D}^{\alpha}X + \frac{\lambda + \mu}{m\_{2}}\mathbf{D}^{\alpha\_{2}}Y - \frac{\kappa\_{1}}{m\_{2}}X + \frac{\kappa\_{1} + \kappa\_{2}}{m\_{2}}Y &= g(\mathbf{x}), \ \mathbf{1} < \beta \le 2, \ \mathbf{0} < \alpha\_{1}, \alpha\_{2} \le 1, \end{aligned} \tag{60}$$

subject to the boundary conditions

$$\begin{aligned} X(\mathbf{0}) &= \mathbf{0}, & X(\mathbf{1}) &= \mathbf{0} \\ Y(\mathbf{0}) &= \mathbf{1}, & Y(\mathbf{1}) &= -\mathbf{1}. \end{aligned} \tag{61}$$

Where the source terms are given by

$$\begin{aligned} f(\mathbf{x}) &= \cos\left(\pi\mathbf{x}\right) - \sin\left(\pi\mathbf{x}\right) - \pi\cos\left(\pi\mathbf{x}\right) - \pi\sin\left(\pi\mathbf{x}\right) - \pi^2\sin\left(\pi\mathbf{x}\right) \\ g(\mathbf{x}) &= \sin\left(\pi\mathbf{x}\right) - 2\cos\left(\pi\mathbf{x}\right) + \pi\cos\left(\pi\mathbf{x}\right) + 2\pi\sin\left(\pi\mathbf{x}\right) - \pi^2\cos\left(\pi\mathbf{x}\right) .\end{aligned} \tag{62}$$

*Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order DOI: http://dx.doi.org/10.5772/intechopen.102054*

Where *λ*, *μ* are damping parameters and *κ*1, *κ*<sup>2</sup> are spring constants and *m*1, *m*<sup>2</sup> are masses of objects and springs are hanged from both ends. We solve this problem under the following values:

$$\begin{aligned} \lambda &= 2, \mu = 2, \kappa\_1 = 2, \kappa\_2 = 2, m\_1 = m\_2 = 2, \\ \text{then calculating the coefficient of the problems as} \\ A\_1 &= \mathbf{1}, A\_2 = -\mathbf{1}, B\_1 = -\mathbf{1}, B\_2 = 2, C\_1 = \mathbf{1}, C\_2 = -\mathbf{1}, D\_1 = -\mathbf{1}, D\_2 = \mathbf{1}. \end{aligned}$$

Let the exact solutions at *α* ¼ *β* ¼ 2, *ν<sup>i</sup>* ¼ *ω<sup>i</sup>* ¼ 1ð Þ *i* ¼ 1, 2 are given by

$$X(\mathfrak{x}) = \sin \left(\pi \mathfrak{x}\right), \qquad Y(\mathfrak{x}) = \cos \left(\pi \mathfrak{x}\right) \tag{63}$$

The above model is approximated for the solutions by our proposed methods (**Figures 6** and **7**). We observed that our method provides best approximate solutions to the problems for small-scale level *M* ¼ 5*:* We also find numerical solutions for fractional values of *α*, *β*, *νi*,*ωi*ð Þ *i* ¼ 1, 2 . We observe that as these orders tend to

**Figure 6.**

*Subplot (a) represents comparison of approximate and exact X x*ð Þ *at M* ¼ 5*. The subplot (b) represents comparison of approximate and exact Y x*ð Þ at *M* ¼ 5 *for example, (4). Blue curve represents exact solution while the red curve represents approximate solution for a* ¼ *b* ¼ 0.

#### **Figure 7***.*

*Subplot (c) represents absolute error* ∣*Xapp* � *Xexact*∣*, and subplot (d) represents their absolute error* ∣*Yapp* � *Yexact*∣ *at M* ¼ 6 *for example, (4) for a* ¼ *b* ¼ 0*.*


**Table 7.** *Maximum absolute error at various values of a*, *b and M for example, (4).*


**Table 8.**

*Maximum absolute error at a* ¼ *b* ¼ 0 *and for various α*, *β*, *ν*1, *ν*2,*ω*1,*ω*<sup>2</sup> *for example, (4).*

their integer values, the solutions are tending to the exact, which prove the applicability of our method in the following tables.

Absolute error for various *a*, *b*, *M* and *α* ¼ *β* ¼ 2, *ν<sup>i</sup>* ¼ *ω<sup>i</sup>* ¼ 1ð Þ *i* ¼ 1, 2 is given in **Table 7**.

We observe from **Table 7** that values of *a*, *b* play important role in the approximate solution. By giving integral value to *a*, *b*, methods give best approximate solution for the given problem. Similarly for *a* ¼ *b* ¼ 0 and for various choices of *α*, *β*, *ν*1, *ν*2,*ω*1,*ω*2, the absolute errors are given in **Table 8**.

From **Table 8**, it is obvious that when the orders of the derivatives approach to their integral values, the error is reducing and the approximate solutions are converging to the exact solutions. This phenomenon indicates that the proposed method is highly accurate.

#### **6. Conclusion and discussion**

In this article, we have studied shifted Jacobi polynomials. Based on these polynomials, we recalled some already existing matrices of fractional order derivative and integrations from the literature as well as we constructed a new operational matrix corresponding to boundary conditions. By means of these operational matrices, we converted the system of fractional differential equations to simple and easily soluble system of algebraic equations. There is no need of Tau-collocation method. The simple algebraic equations were easily solved, and the result were plotted. From the plot, we observe that our proposed method is highly accurate and can be applied to a variety of problems of fractional order ordinary as well as partial differential equations. We also compared our results to the exact solutions and observed that our method gave satisfactory results. The proposed method can easily *Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order DOI: http://dx.doi.org/10.5772/intechopen.102054*

and accurately can be applied to a variety of problems of applied mathematics, engineering, and physics, etc.

### **Acknowledgements**

The authors Aziz Khan, Kamal Shah and Thabet Abdeljawad would like to thank Prince Sultan university for the support through TAS research lab.

### **Mathematics Subject Classification:**

26A33; 34A08; 34K37.

### **Author details**

Kamal Shah1,2\*,Eiman2 , Hammad Khalil<sup>3</sup> , Rahmat Ali Khan<sup>2</sup> and Thabet Abdeljawad1,4

1 Department of Mathematics, University of Malakand, Khyber Pakhtunkhwa, Pakistan

2 Department of Mathematics and Sciences, Prince Sultan University, Riyadh, Saudi Arabia

3 Department of Mathematics, University of Education Lahore, Attock, Punjab, Pakistan

4 Department of Medical Research, China Medical University, Taichung, Taiwan

\*Address all correspondence to: kamalshah408@gmail.com

© 2022 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/ by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

### **References**

[1] Podlubny I. Fractional Differential Equations, Mathematics in Science and Engineering. New York: Academic Press; 1999

[2] Lakshmikantham V, Leela S, Vasundhara J. Theory of Fractional Dynamic Systems. Cambridge, UK: Cambridge Academic Publishers; 2009

[3] Hilfer R. Applications of Fractional Calculus in Physics. Singapore: World Scientific; 2000

[4] Rossikhin YA, Shitikova MV. Applications of fractional calculus to dynamic problems of linear and nonlinear hereditary mechanics of solids. Applied Mechanics Reviews. 1997;**50**:15-67

[5] Din A, Li Y, Khan FM, et al. On Analysis of fractional order mathematical model of Hepatitis B using Atangana-Baleanu Caputo (ABC) derivative. Fractals. 2021;**30**:2240017

[6] Lakshmikantham V, Leela S. Naguma-type uniqueness result for fractional differential equations. Nonlinear Anal. 2009;**71**:2886-2889

[7] Miller KS, Ross B. An Introduction to the Fractional Calculus and Fractional Differential Equations. New York: Wiley; 1993

[8] Goodrich C. Existence of a positive solution to a class of fractional differential equations. Computers & Mathematcs with Applications. 2010;**59**: 3889-3999

[9] Din A, Li Y. The extinction and persistence of a stochastic model of drinking alcohol. Results in Physics. 2021;**28**:104649

[10] Ahmad B, Nieto JJ. Existence result for a coupled system of non linear

fractional differential equations with three point boundary conditions. Computers & Mathematcs with Applications. 2009;**58**:1838-1843

[11] Din A, Li Y. Lévy noise impact on a stochastic hepatitis B epidemic model under real statistical data and its fractalfractional Atangana-Baleanu order model. Physica Scripta. 2021;**96**(12): 124008

[12] Ray SS, Bera RK. Solution of an extraordinary differential equation by Adomian decomposition method. Journal of Applied Mathematics. 2004; **4**:331-338

[13] Hashim I, Abdulaziz O, Momani S. Homotopy analysis method for fractional IVPs. Commun. Nonl. Sci. Numer. Simul. 2009;**14**:674-684

[14] Din A, Li Y. Stationary distribution extinction and optimal control for the stochastic hepatitis B epidemic model with partial immunity. Physica Scripta. 2021;**96**(7):074005

[15] Khalil H, Khan RA. A new method based on Legendre polynomials for solutions of the fractional twodimensional heat conduction equation. Computers & Mathematics with Applications. 2014;**67**:1938-1953. DOI: 10.1016/j.camwa.2014.03.008

[16] Doha EH, Bhrawy AH, Ezz-Eldien SS. Efficient Chebyshev spectral methods for solving multi-term fractional orders differential equations. Applied Mathematical Modelling. 2011; **35**(12):5662-5672

[17] Doha EH, Bhrawy AH, Ezz-Eldien SS. A new Jacobi operational matrix: An application for solving fractional differential equations. Applied Mathematical Modelling. 2012;**36**: 4931-4943

*Using Shifted Jacobi Polynomials to Handle Boundary Value Problems of Fractional Order DOI: http://dx.doi.org/10.5772/intechopen.102054*

[18] Saadatmandi A, Dehghan M. A new operational matrix for solving fractional-order differential equations. Computers & Mathematcs with Applications. 2010;**59**(3):1326-1336

[19] Singh AK, Singh VK, Singh OP. The Bernstein operational matrix of integration. Appl. Math. Sciences. 2009; **3**(49):2427-2436

[20] Rehman M u, Khan RA. the legender wavelet method for solving fractional differential equation. Commun. Nonlin. Sci. Numer. Simul. 2011;**10**:4163-4173. DOI: 10.1016/ j.cnsns.2011.01.014

[21] Doha EH, Bhrawy AH, Saker MA. On the Derivatives of Bernstein Polynomials, An application for the solution of higher even-order differential equations. Boun. Val. Prob. 2011;**16**:829543. DOI: 10.1186/1687– 2770-2011-16

[22] Khalil H, Khan RA. A new method based on legender polynomials for solution of system of fractional order partial differential equation. International Journal of Computer Mathematics. 2014;**67**:1938-1953. DOI: 10.1080/00207160.2014.880781

[23] Tan B, Boyd JP. Stability and long time evolution of the periodic solutions of the two coupled nonlinear Schrödinger equations. Chaos. Soli. Frac. 2001;**12**:721-734

[24] Radhakrishnan R, Sahadevan R, Lakshmanan M. Integrability and singularity structure of coupled nonlinear Schrödinger equations. Chao. Soli. Fract. 1995;**5**:2315-2327

[25] Bakkyaraj T, Sahadevan R. On solutions of two coupled fractional time derivative Hirota equations. Nonlinear Dynamics. 2014;**77**:1309-1322

[26] Rienstra SW. Non-Linear free vibrations of coupled spans of overhead transmission lines. In: Manley J et al., editors. Proceedings of the Third European Conference on Mathematics in Industry, August 27–31 1990 Glasgow. Alphen aan den Rijn, Netherlands: Kluwer Academic Publishers and B.G. Teubner Stuttgart; 1990. pp. 133-144

[27] Din A et al. Caputo type fractional operator applied to Hepatitis B system. Fractals. 2021;**31**:2240023

[28] Bhrawy AH, Alofi AS. A Jacobi-Gauss collocation method for solving nonlinear Lane-Emden type equations. Commun Nonlinear. Sci. Numer. Simulat. 2012;**17**:62-70. DOI: 10.1016/ j.cnsns.2011.04.025

[29] Kilbas AA, Srivastava H, Trujillo J. Theory and application of fractional differential equations, North Holland Mathematics Studies. Vol. 204. Amsterdam: Elseveir; 2006

[30] Doha EH, Bhrawy AH, Hafez RM. On shifted Jacobi spectral method for high-order multi-point boundary value problems. Commun. Nonl. Sci. Num. Simul. 2012;**17**(10):3802-3810

[31] Bhrawy AH, Alghamdi MA. Numerical solutions of odd order linear and nonlinear initial value problems using shifted Jacobi spectral approximations. Abst. Appl. Anal. 2012; **2012**:364-360

## *Edited by Kamal Shah*

This book provides a broad overview of recent developments in polynomials and their applications. It includes eight chapters that address such topics as characteristic functions of polynomials, permutations, Gončarov polynomials, irreducible factors, polynomial regression algorithms, and the use of polynomials in fractional calculus, and much more.

Published in London, UK © 2022 IntechOpen © WhataWin / iStock

Recent Advances in Polynomials

Recent Advances in

Polynomials

*Edited by Kamal Shah*