Abstract
In this paper, a new approach for solving the system of fractional integro-differential equation with weakly singular kernels is introduced. The method is based on a class of symmetric orthogonal polynomials called shifted sixth-kind Chebyshev polynomials. First, the operational matrices are constructed, and after that, the method is described. This method reduces a system of weakly singular fractional integro-differential equations (WSFIDEs) by the collocation method into a system of algebraic equations. Thereupon, an upper error bound for the proposed method is determined. Finally, some numerical examples are prepared to test the accuracy and efficiency of the presented method.
1. Introduction
The study of fractional calculus has applications and popularity in various and wide fields of biology, physics, and fluid mechanics. Fractional calculus is actually integration and differentiation of arbitrary orders [1–4]. In various problems of physics and engineering, the fractional differential equations have been proved to be valuable tools in modeling of many phenomena [5, 6]. As we know, many mathematical models of real phenomena (arising in engineering and physics) are described as linear or nonlinear systems. It is worth mentioning that with the development of fractional calculus, the behavior of many systems can be described using the fractional differential and fractional integro-differential system [7, 8]. In recent years, systems of the fractional differential and integral equations are the subjects of extensive study due to their frequent appearance in many engineering and scientific disciplines [9–11]. However, most of the fractional-order equations and integral equations do not have analytic solutions or are hard to find. So, it is essential to find numerical methods to get approximate or exact solutions of a system of integro-differential equations. So far, researchers have utilized diverse numerical methods for solving a system of fractional integro-differential equations. In [12], the homotopy perturbation method was proposed for solving linear and nonlinear systems of fractional integro-differential equations. Heydari et al. have used the Chebyshev wavelet method for solving a class of systems of nonlinear fractional singular Volterra integro-differential equation in [13]. In the next year, for the first time, the hybrid functions composed of the Block-pulse functions and Bernoulli polynomials were applied for problems with fractional-order differential equations in [14]. Also, a novel technique based on iterative refinement was presented to analytically approximate a system of linear fractional integro-differential equations [15]. In 2018, Hesameddini and Shahbazi developed the concept of [14] and used it to solve the FDIE system in [16]. Also, Xie and Yi presented the simple and fast method based on the Block-pulse function to solve a nonlinear system of fractional Volterra-Fredholm integro-differential equations in the same year [17]. In [18], the authors implemented the new Jacobi operational matrices to reduce the complexity of calculations to solve WSFDIEs. Next, the Haar wavelet method was employed to solve a coupled system of FIDEs, and the Muntz-Legendre wavelets were introduced to solve FIDVFEs in [19, 20]. Also, the other authors applied the Chebyshev Pseudo spectral method for solving fractional-order nonlinear system of Volterra integro-differential equations and a least square collocation Chebyshev technique for solving a system of linear fractional integro-differential equations [21, 22]. In this paper, we consider the following system of weakly singular integro-differential equations:
where are the unknown functions, , , and are continuous operators and functions that satisfy Lipschitz conditions, and is the Caputo fractional derivative operator where . The parameters such that , , and . Moreover, are known and sufficiently smooth functions.
As usual, a way for solving functional equations is to express the solution as a linear combination of the so-called basis functions. In most researches, various polynomials such as the Legendre, Chebyshev, Taylor, Hermit, and Bernstein are used as basis functions. Among all of them, the Chebyshev polynomials are the most important in the analysis and numerical analysis. Chebyshev polynomials are orthogonal on the interval and have good properties that are used widely in the approximation of the functions. For this reason, many studies are done based on the different kinds of Chebyshev polynomials. In [23], Masjed-Jamei introduced two half-trigonometric orthogonal Chebyshev polynomials, and he named them as the Chebyshev polynomials of the fifth and sixth kinds. The basic formulas and properties of this class of polynomials are displayed in [24, 25]. Up to now, many researchers have used various kinds of Chebyshev polynomials for the fractional-order differential and integro-differential equations (see [26, 27]). However, there are only a few works that have used the sixth-kind Chebyshev polynomials. The main aim of this work is to introduce these polynomials as a new basis function for solving WSFDIEs. In the current paper, we apply the orthogonal shifted sixth-kind Chebyshev polynomials together with the collocation method for solving a system of weakly singular integro-differential equations with fractional derivatives that, to the best of our knowledge, is proposed here for the first time. For solving these equations, we derive the fractional operational matrices of fractional and integer orders and the product operational matrix, as well. Also, we introduce an operational matrix to approximate the integral term that has the weakly singular kernel in Equation (1). As far as we can tell, this operational matrix is presented for the first time. By substituting appropriate approximations in Equation (1), the original equations convert into algebraic equations that each of the equations of algebraic systems is collocated at roots of the th shifted sixth-kind Chebyshev polynomials (SSKCPs). By solving these algebraic systems, the approximate solutions of the original system are obtained. Although the calculation of the operational matrices may be complicated, we show that the obtained results are equal to other methods or are even more accurate. Implementing these matrices leads to a decrease in the number of required computations, and therefore, the computation time will be reduced.
The rest of the paper is organized as follows. In section 2, some essential preliminaries are mentioned briefly. Section 3 is devoted to constructing the operational matrices of SSKCPs. The proposed numerical procedure is described in Section 4. The error analysis of the proposed method is discussed in Section 5. Some numerical applications are indicated in Section 6, and conclusions are presented in Section 7.
2. Preliminaries and Notations
In this section, we recall some definitions and properties of fractional integral and derivative operators which will be used later. After that, some necessary definitions and fundamental properties of the shifted sixth-kind Chebyshev polynomials are reviewed briefly.
2.1. Some Essentials of the Fractional Calculus
Definition 1. The Riemann-Liouville fractional integral operator of order is given by [2]
Definition 2. Let , , , and be a real-valued continuous function defined on . Then, the Caputo fractional derivative of order is defined by [2] where is the Gamma function as
The last integral is often called the Beta integral. For , the Riemann-Liouville integral and Caputo fractional derivative operators satisfy the following properties: (1)(2)(3)(4)(5)
2.2. Shifted Sixth-Kind Chebyshev Polynomials and Their Properties (SSKCPs)
Definition 3. The sixth-kind Chebyshev polynomials are orthogonal functions on the interval and can be determined with the following recursive formula [23, 24]:
Definition 4. The shifted sixth-kind Chebyshev polynomials on is defined by [23, 24] These polynomials have the following explicit analytic form: where
Moreover, the shifted polynomials are orthogonal on with respect to the weight function in the sense that
Now, let ; then, can be approximated in terms of as where where the coefficients are given by and is defined in Equation (10). Similarly, any continuous two-variable function, , defined on can be approximated by means of the double-shifted sixth-kind Chebyshev polynomials as where is a matrix, and its entries are given by
3. Operational Matrices of SSKCPs
In this section, the formulas of operational matrices with the fractional order will be derived for the sixth-kind Chebyshev polynomials. The following are the required lemmas.
Lemma 5. If , , then we have
Proof. From the properties of the orthogonal polynomials, if , we have Hence, we suppose . The lemma can be easily proved by the integration of the analytic form of SSKCPs in Equation (7).
Theorem 6. Let be the SSKCP vector given by Equation (12) and ; then, where is the operational matrix of the fractional integration of the order in the Riemann-Liouville sense which is defined by are given by
Proof. By applying the Riemann-Liouville integral operator to the SSKCPs’ analytic form, we have Now, we can express in terms of the shifted sixth-kind Chebyshev polynomials as follows: where the coefficients are given by According to Lemma 5, we can rewrite Equation (22) as where is given in Equation (20). Rewriting the last relation in the vector form gives This leads to the desired result.
In the following, some useful and applicable lemmas are presented to get the Chebyshev operational matrix of product.
Lemma 7. If and are th and th shifted sixth-kind Chebyshev polynomials, respectively, then we can write the product of and as
Proof. See [18].
Lemma 8. If , , and are th, th, and th shifted sixth-kind Chebyshev polynomials, then where is obtained by Lemma 7.
Proof. According to Lemma 7, we can write Then, The value of the last integral is obtained by Lemma 5.
Assuming that is a vector, we have where is a matrix called the product operational matrix. The next theorem presents a general form for entries of the matrix .
Theorem 9. The entries of the matrix in Equation (31) are as follows: where is obtained by Lemma 8 and is the element of the vector .
Proof. See [18].
In the following, we get an approximation for the integral part with the singular kernel in Equation (1). Before that, we present a theorem.
Theorem 10. The following relation is determined for :
Proof. By performing Equation (33) and the substitution of into Equation (33) and then using the definition of the Beta function, we obtain
Now, we present an approximation for the integral part with a weakly singular kernel. For this purpose, see the following theorem.
Theorem 11. Suppose that is a continuous function on the interval and and where and are defined by Equation (12), then we have where is a matrix as follows: and its entries are determined as follows: where the quantities and are introduced by relation (8) and Lemma 5.
Proof. By the definition of vector and Lemma 5, we can write By applying Theorem 10, we have Now, we approximate in terms of SSKCPs as follows: where is obtained by applying Lemma 5. Thus, Equation (39) is then represented as follows:
4. Numerical Procedure
We consider the system of fractional integro-differential equation with weakly singular kernels described in Equation (1). To solve the system, we approximate the function in a matrix form
Then, according to the initial conditions of the problem, we can approximate the known function as follows:
Using Equations (43) and (44), we compute an approximate for : where is the integral operational matrix presented in Theorem 6. We have the following approximations for the rest of the system:
Using Theorems 9 and 11, we have
Now, by substituting Equations (43)-(47) into Equation (1), we obtain
Each equation of algebraic system (48) is collocated at roots of the shifted sixth-kind Chebyshev polynomials. Thus, an algebraic system, including equations, is acquired. By solving the resultant algebraic system, we can obtain an approximation for unknown vectors , and by substituting the vector into Equation (45), we obtain an approximation for .
5. Error Analysis
In this section, we prove some theorems. Then, we obtain an upper error bound for the approximation error. For this aim, we need the following norms: where is a square integrable function on the interval and is a vector.
Theorem 12. Suppose that is an approximation in SSKCPs to the continuous function on the interval . Then, the coefficients , for , are bounded as where denotes the maximum value of on the interval .
Proof. Using Equations (7) and (9) for we have Since is a continuous function on the interval , so it is bounded and there is a constant such that Using Equations (51) and (52), inequality (50) is deduced.
Theorem 13. Suppose that is a continuous function and is an approximation to in terms of SSKCPs. Then, the error bound can be achieved as follows: where
Proof. Assume is an arbitrary function. So, and have the following forms using SSKCP series: so, Using Equations (9) and (56) and Theorem 12, we have
Theorem 14. Suppose that the continuous two-variable function is approximated on the interval in terms of SSKCPs as ; then, the coefficients , for , can be bounded as follows: where denotes the maximum value of on the interval .
Proof. Using Equations (7) and (15), we have Since is a bounded and continuous function on the interval , so there is a constant such that Using Equations (59) and (60), Theorem 14 is proved.
Theorem 15. Suppose that is a continuous function with two variables and is the approximation to using SSKCPs. Then, the error bound can be obtained as where
Proof. Suppose that is an arbitrary function. SSKCP series of and its approximation in terms of SSKCPs have the following form: Thus, Using Equations (9) and (64) and Theorem 14, we conclude
In the following theorem, we obtain an upper error bound of the proposed method. First, suppose that for the there exist positive constants such that the following Lipschitz conditions hold. (1)(2)
Theorem 16. Suppose that is a set of approximate solutions obtained from the SSKCP collocation method, is the set of exact solutions of system (1), is the error vector of approximate solutions, and also denote the residual functions associated to the approximate solutions that are named perturbation terms. Assume that Hypotheses (1) and (2) are satisfied; then, a bound for the method error can be achieved as
Proof. First, we apply the Riemann-Liouville integral operator on Equation (1) and obtain the following equation: where We can write the approximate equation of Equation (67) as follows: where is the perturbation term. We subtract Equation (69) from Equation (67) and obtain the following result: First, we obtain a bound for the perturbation term, so by taking the -norm on Equation (70), we get the following inequality: Using Hypothesis 1 and Theorem 13, we have where Noting that and are continuous and known functions, thus, there are constants and such that From Hypothesis 2, Theorem 13, Theorem 15, and Equation (74), the following inequality is obtained: where Here, and can be introduced as below: