A New Operational Matrix of Fractional Derivatives to Solve Systems of Fractional Differential Equations via Legendre Wavelets

This paper introduces a new numerical approach to solving a system of fractional differential equations (FDEs) using the Legendre wavelet operational matrix method (LWOMM). We first formulated the operational matrix of fractional derivatives in some special conditions using some notable characteristics of Legendre wavelets and shifted Legendre polynomials. Then, the system of fractional differential equations was transformed into a system of algebraic equations by using these operational matrices. At the end of this paper, several examples are presented to illustrate the effectivity and correctness of the proposed approach. Comparing the methodology with several recognized methods demonstrates that the advantages of the Legendre wavelet operational matrix method are its accuracy and the understandability of the calculations.


Introduction
Differential and integral operators are the basis of mathematical models, and they are also used as a means of understanding the working principles of natural and artificial systems.Therefore, differential and integral equations are of great importance both theoretically and practically.Such equations have a wide range of applications, including in the physical sciences (such as in physics and engineering) as well as in social science.Systems of differential equations, as differential equations, are often used in issues such as theories of elasticity, dynamics, fluid mechanics, oscillation, and quantum dynamics.
Interest in differential and integral operators has led to the exploration of fractional differential and integral operators by examining these issues further in depth.Owing to a question, the origin of fractional calculus arose in a message from Leibniz to L'Hôpital in 1695.Fractional calculus has received attention in recent years due to its ability to simplify numerous physical, engineering, and economics phenomena such as the fluid dynamic traffic model, damping laws, continuum and statistical mechanics, diffusion processes, solid mechanics, control theory, colored noise, viscoelasticity, electrochemistry, and electromagnetism, among others.
Because a variety of solutions of fractional differential equations (FDEs) cannot be found analytically, numerical and approximate methods are needed.There are a lot of techniques that have been studied by many researchers in solving FDEs and the system of such equations numerically.Some of these techniques are the Adomian decomposition method presented by Song and Wang [1], the collocation method, the improved operational matrix method [2][3][4], the perturbation iteration method introduced by Şenol and Dolapçı [5], the computational matrix method illustrated by Khader et al. [6], the differential transform method demonstrated by Ertürk and Momani [7], the variational iteration method, the Laplace transform method given by Gupta et al. [8], and the fractional complex transform method studied by Ghazanfari and Ghazanfari [9], among others.Kilbas et al. [10] inclusively examined fractional differential and fractional integro-differential equations.In addition, numerical solutions of FDEs and the system of such equations have been presented using the Legendre polynomial operational matrix method [11], Bernstein operational matrix method [12], Genocchi operational matrix method [13], Jacobi operational matrix method [14], Chebyshev wavelet operational matrix method [15], polynomial least squares method (PLSM) [16], Legendre wavelet-like operational matrix method (LWPT) [17], and the Genocchi wavelet-like operational matrix method [18].
This paper focuses on the numerical analysis of a system of fractional order differential equations using the Legendre wavelet operational matrix method.The most important advantage of the proposed method is that it presents an understandable procedure to reduce FDEs and the system of such equations to a system of algebraic equations.First, we begin by presenting some basic definitions and fundamental relations in Sections 2 and 3, respectively.Then, in Section 4, the operational matrix of the fractional derivate is natively formulated to linear and nonlinear systems of fractional differential equations.Section 5 presents five illustrative examples that were tested with the introduced method.Finally, the last section includes the conclusions.

Basic Definitions
The Liouville_Caputo fractional_order derivative, shifted Legendre polynomials, and Legendre wavelets are defined below [19,20].Definition 1.The Liouville_Caputo fractional derivative of u is defined as [19] Some characteristics of the Liouville_Caputo fractional derivative are as follows: where C is a constant.In addition, there is in which α and α respectively imply that the largest integer is less than or equal to α, and the smallest integer is greater than or equal to α.The Liouville_Caputo fractional order derivative is a linear operation of the integer order derivative where η and ζ are constant.
Definition 2. Let a and b respectively be the parameters of dilation and translation of a single function called the mother wavelet.If a and b change continuously, then we obtain the following family of continuous wavelets [21,22]: Definition 3. Let P m (t) imply the shifted Legendre polynomials of order m.Then P m (t) can be formulated as [21] P and the orthogonality condition is Definition 4. Let n and k be any positive integer, m be the order of shifted Legendre polynomials, and t be the normalized time.Then the Legendre wavelets ψ nm (t) = ψ(k, n, m, t) are defined on the interval [0, 1]  by [21,22].

Fundamental Relations
Saadatmandi and Dehghan [11] derived the operational matrix of a fractional derivative by using shifted Legendre polynomials.In this section, we show how we derived the Legendre wavelet operational matrix of fractional derivatives in some special conditions by drawing from Saadatmandi and Dehghan [11].Additionally, the theorem and corollary related to the Legendre wavelet operational matrix of derivatives illustrated by Mohammadi [21] are cited here as follows.
Theorem 1.Let ψ(t) be the Legendre wavelet vector introduced in Equation ( 8).Then ψ(t) is expressed as [21,22] dψ where D is the 2 k (M + 1) operational matrix of the derivative, which can be stated as where U is an (M + 1)(M + 1) matrix and its (r, s)th element is written as Corollary 1.Using Equation ( 12), the operational matrix for the nth derivative can be stated as [21] d n ψ(t) where D n is the nth power of matrix D.
Theorem 2. Let ψ(t) be the Legendre wavelets vector introduced in Equation ( 8).Supposing that k = 0 and α > 0, then where D (α) is the (M + 1)x(M + 1) operational matrix of the fractional derivative of the order α > 0, N − 1 < α ≤ N in the Liouville_Caputo sense and can be stated as where ξ r,s,h is written as Consider in D (α) that the first α rows are all zero.

Solving Systems of Fractional Order Differential Equations
In this section, the Legendre wavelet operational matrix method was implemented to obtain the numerical solution of the system of fractional order differential equations.Consider a system of fractional differential equations as follows: where U i is a linear/nonlinear function of t, u 1 , u 2 , . . ., u m , D η i is the derivative of u i with the order of η i in the Liouville-Caputo sense and N − 1 ≤ η i < N, and they are subjected to the initial conditions First of all, approximating u 1 (t), u 2 (t), . . . ,u m (t) and D η 1 u 1 (t), D η 2 u 2 (t), . . . ,D η n u m (t), we obtain . . .
where C i , i = 1, 2, . . ., m is an unknown vector and ψ(t) is the vector introduced in Equation ( 8).If we utilize Equation ( 17), then we have Substituting Equations ( 29) and (30) into Equation (27), we obtain If U i is a linear function of t, u 1 , u 2 , . . ., u m , then we produce 2 k (M + 1) − mn linear equations by implementing Also, by substituting the initial conditions in Equation (28) into Equation (30), then we obtain . (33) A 2 k (M + 1) set of linear equations is generated by combining Equations ( 32) and (33).The solution of these linear equations can be obtained for unknown coefficients of the vector C. Consequently, u 1 (t), u 2 (t), . . . ,u m (t), introduced in Equation ( 27), can be computed.

Illustrative Examples
In this section, to show the applicability and powerfulness of the introduced method, we present the solutions to five linear and nonlinear systems of fractional order differential equations.
Example 1.We first considered the following linear system of fractional differential equations [7,8]: The exact solution of this system when α = 1 is known to be u(t) = e t sin t, v(t) = e t cos t.
This example was examined for M = 2, k = 0, and α = 0.9, 0.7, 0.5.When the obtained results were matched against the exact solution when α = 1, as demonstrated in Figure 1, we can clearly observe that when α approached 1, our results approached the exact solution.We also solved this problem by using Legendre polynomial operational matrix method (LPOMM), and we compared the results with the LWOMM.The numerical computations for u(t) and v(t) when α = 0.9 are revealed in Tables 1 and 2.    Example 2. We considered the following nonlinear system of fractional differential equations [13]: The exact solution of this system is known to be ( ) ( ) Using the parameters 3 M = and 0 k = , we applied both the proposed method and the LPOMM to solve this problem and show that our approach is more efficient and useful.Our numerical results supported the idea that our solution approaches the exact solution more than the approximate solution LPOMM.Comparisons of the approximate and exact solutions are presented in Tables 3 and 4. Example 2. We considered the following nonlinear system of fractional differential equations [13]: The exact solution of this system is known to be Using the parameters M = 3 and k = 0, we applied both the proposed method and the LPOMM to solve this problem and show that our approach is more efficient and useful.Our numerical results supported the idea that our solution approaches the exact solution more than the approximate solution LPOMM.Comparisons of the approximate and exact solutions are presented in Tables 3 and 4. Example 3. We considered the following nonlinear system of fractional differential equations with the initial conditions The exact solution of this system when α = 1 is known to be The parameters M = 2, k = 0, and α = 0.5, 0.7, 0.9 were utilized.A comparison of our results and the exact solution when α = 1 is displayed in Figure 2. The figures support that when α approximated 1, our results approximated the exact solution.We also solved this problem by using the LPOMM, and compared the results to the LWOMM.Finally, we present the numerical computations for u(t) and v(t) when α = 0.9 in Tables 5 and 6.
The exact solution of this system when α = 1 is known to be This example was analyzed for M = 3, k = 0, and α = 0.9, 0.7, 0.5.When the obtained results were matched against the exact solution when α = 1, as demonstrated in Figure 3, we can clearly observe that when α approached 1, our results approached the exact solution.We also solved this problem by using the LPOMM, and compared the results with the LWOMM.The numerical computations for u(t) and v(t) when α = 0.99 are also revealed in Tables 7 and 8.

Conclusions
In this paper, a system of fractional_order differential equations was examined by drawing from a new operational matrix of the fractional derivative in some special conditions.We also systematized a very operational algorithm in order to attain the solution of the linear and nonlinear systems of fractional differential equations in Maple.All numerical results and graphical presentations generated by Maple affirmed that the Legendre wavelet operational matrix method is very effective and applicable.As the next step, the method introduced in this paper can be applied to fractional partial differential equations and the system of such equations, fractional integral equations and the system of such equations, and fractional integro-differential equations.These equations are at least as important as fractional differential equations and they are very significant in science, engineering, and technology.

Figure 1 .
Figure 1.Comparison of our solutions and the exact solution when α = 0.9, 0.7, 0.5 in Example 1: (a) Our solution u(t); and (b) Our solution v(t).

Figure 2 .
Figure 2. Comparison of our solutions to the exact solution when α = 0.9, 0.7, 0.5 for Example 3: (a) Our solution u(t); and (b) Our solution v(t).

Figure 4 .
Figure 4. Comparison of our solutions to the approximate solution LWPT and the approximate solution PLSM when α = 0.98 for Example 5: (a) Our solution u(t); and (b) Our solution v(t).

Table 1 .
Numerical solutions of u(t) when α = 0.9 attained by the introduced method and the LPOMM for Example 1.

Table 2 .
Numerical solutions of v(t) when α = 0.9 attained by the introduced method and the LPOMM for Example 1.

Table 2 .
Numerical solutions of v(t) when

Table 3 .
The numerical results attained by using the introduced method in comparison to the approximate solution LPOMM and the exact solution u(t) in Example 2.

Table 4 .
The numerical results attained by using the introduced method in comparison with the approximate solution LPOMM and exact solution v(t) in Example 2.

Table 5 .
Our solutions u(t) when α = 0.9 attained by the presented method and the LPOMM for Example 3.

Table 6 .
Our solutions v(t) when α = 0.9 attained by the presented method and the LPOMM for Example 3.Example 4. We considered the following nonlinear system of FDEs with initial conditions[13]

Table 10 .
Numerical solutions of v(t) when α = 0.98 obtained by the introduced method, the LWPT, and the PLSM for Example 5.

Table 10 .
Numerical solutions of v(t) when obtained by the introduced method, the LWPT, and the PLSM for Example 5.