A Novel Iterative Method for Solving Systems of Fractional Differential Equations

Abstract

The Reconstruction of Variational Iteration Method (RVIM) technique has been successfully applied to obtain solutions for systems of nonlinear fractional differential equations: 𝐷^𝛼_∗𝑥_𝑖(𝑡)=𝑁_𝑖(𝑡,𝑥₁,…,𝑥_𝑛), 𝑥_𝑖^(𝑘)=𝑐^𝑖_𝑘,  0≤𝑘≤[𝛼_𝑖], 1≤𝑖≤𝑛, where 𝐷^𝛼_∗ denote Caputo fractional derivative. The RVIM, for differential equations of integer order is extended to derive approximate analytical solutions for systems of fractional differential equations. Advantage of the RVIM, is simplicity of the computations and convergent successive approximations without any restrictive assumptions or transform functions. Some illustrative examples are given to show the validity of this method for solving linear and nonlinear systems of fractional differential equations.

1. Introduction

In recent years, the fractional differential equations have received remarkable attention. Differential equations of fractional order have been found to be effective to describe some physical phenomena such as rheology, fluid flow, diffusive transport, electrical network, and electromagnetic theory [1–4]. There are different methods to solve the fractional differential equations. Some of the recent analytic methods for solving a system of nonlinear fractional differential equations are the Adomian decomposition method (ADM) [5–8], differential transform method [9], and Variational Iteration method (VIM) [10].

The differential transform method was first applied in engineering in 1986 [11]. Ertürk and Momani introduced a new application of the differential transform method to provide an approximate solution for systems of fractional differential equations [9]. For this propose He developed the Variational Iteration Method (VIM) in 1999 [10]. In this method, the solution is approximated at first iteration by using the initial conditions. A correction functional is established by the general Lagrange multiplier which can be identified optimally via the variational theory. Although a number of useful attempts have been made to solve fractional equations via the VIM, the problem has not yet been completely resolved; that is, most of the previous work avoid the term of fractional derivative and handle them as restricted variations and they cannot identify the fractional Lagrange multipliers explicitly in the correction function. Hesameddini and Latifizadeh proposed a new alternative approach to derivation of the Variational Iteration formulations using the Laplace transform for solving linear and nonlinear ordinary differential equations which was called the Reconstruction of Variational Iteration Method [12]. This method does not use a Lagrange multiplier.

Partial differential equations of fractional order are often very complicated to be exactly solved and even if an exact solution is obtainable, the required calculations may be too complicated to be practical, or it might be difficult to interpret the outcome.

In this work, we extend the Reconstruction of Variational Iteration Method to solve systems of fractional differential equations. The aim of this work is to present an alternative approach based on RVIM to find the solution for linear and nonlinear system of fractional differential equations. The efficiency and accuracy of RVIM are demonstrated through several test examples.

2. Preliminaries and Notations

In this section, some necessary definitions and mathematical preliminaries of the fractional calculus theory which are used further in this paper will be presented.

Definition 1. Let 𝐶[𝑎,𝑏] denotes the space of continuous functions defined on [𝑎,𝑏] and 𝐶^𝑛[𝑎,𝑏] denotes a class of all real valued functions defined on [𝑎,𝑏] which have continuous 𝑛th order derivative.

Definition 2. Let 𝑓∈𝐶[𝑎,𝑏] and 𝛼≥0; then the expression 𝐼^𝛼𝑓(𝑡)=1Γ(𝛼)^𝑡₀(𝑡−𝜏)^𝛼−1𝑓(𝜏)𝑑𝜏, (1) is called the Riemann-Liouville integral of order 𝛼.

Definition 3. The fractional derivative of 𝑓(𝑥) in the Caputo sense is defined as 𝐷^𝛼_∗𝑓(𝑡)=𝐼^𝑚−1𝑓^(𝑚)(𝑡)=1Γ()×∫𝑚−𝛼_𝑡0(𝑡−𝜏)^{𝑚−𝛼−1}𝑓^(𝑚)(𝜏)𝑑𝑑𝜏,𝑚−1<𝛼<𝑚,_𝑚𝑑𝑥^𝑚𝑓(𝑡),𝛼=𝑚, (2) for 𝑚∈𝑁 and 𝑓∈𝐶^𝑛[𝑎,𝑏].

Note that 𝐼^𝛼𝐼^𝛽𝑓(𝑡)=𝐼^𝛼+𝛽𝑓(𝑡)𝐼,𝛼,𝛽≥0,_𝛼𝑡^𝛾=Γ𝛾+1Γ𝑡𝛾+𝛼+1^𝛾+𝛼,𝛼>0,𝛾>−1,𝑡>0. (3)

Definition 4. Given a function 𝑓(𝑡) defined for all  𝑡≥0, the Laplace transform of 𝑓 is the function 𝐹 defined as follows: 𝐹(𝑠)=𝑙𝑓(𝑡;𝑠)=^∞₀𝑒^−𝑠𝑡𝑓(𝑡)𝑑𝑡, (4) for all values of 𝑠 for which the improper integral converges.

Definition 5. The function 𝑓 is said to be of exponential order as 𝑡→+∞ if there exist nonnegative constants 𝑀,  𝑐, and 𝑇 such that 𝑓(𝑡)≤𝑀𝑒^𝑐𝑇for𝑡≥𝑇. (5)

Definition 6. The function 𝑓(𝑡) is said to be piecewise continuous on the bounded interval 𝑎≤𝑡≤𝑏 provided that [𝑎,𝑏] can be subdivided into finitely many abutting subintervals in such a way that(1)𝑓 is continuous in the interior of each of these subintervals,(2)𝑓(𝑡) has a finite limit as 𝑡 approaches each endpoint of each subinterval from the interior.

Theorem 7 (existence of the laplace transforms). If the function 𝑓 is piecewise continuous for 𝑡≥0 and is of exponential order as 𝑡→+∞, then its Laplace transform 𝐹(𝑠) exists. More precisely, if 𝑓 is piecewise continuous and satisfies the condition (5), then 𝐹(𝑠) exists for all 𝑠>𝑐.

Definition 8. Let the functions 𝑓(𝑡) and 𝑔(𝑡) be defined for 𝑡≥0; then the convolution of them is denoted by (𝑓∗𝑔)(𝑡) and is defined as the following integral: 𝑓∗𝑔(𝑡)=^𝑡₀𝑓(𝜏)𝑔(𝑡−𝜏)𝑑𝜏. (6) In other words, if 𝑙{𝑓(𝑡)}=𝐹(𝑠),  𝑙{𝑔(𝑡)}=𝐺(𝑠), then 𝑙{(𝑓∗𝑔)(𝑡)}=𝐹(𝑠)𝐺(𝑠). Or equivalently, 𝑙{∫^𝑡₀𝑓(𝜏)𝑔(𝑡−𝜏)𝑑𝜏}=𝐹(𝑠)𝐺(𝑠). Therefore, the inverse Laplace transform will be defined as 𝑙⁻¹{𝐹(𝑠)𝐺(𝑠)}=^𝑡₀𝑓(𝜏)𝑔(𝑡−𝜏)𝑑𝜏. (7)

Definition 9. The Laplace transform of the Caputo fractional derivative 𝐷^𝛼_∗𝑓(𝑡) is given by 𝑙𝐷^𝛼_∗𝑓(𝑡);𝑠=𝑠^𝛼𝐹(𝑠)−^𝑚−1_𝑘=0𝑠^{(𝛼−𝑘−1)}𝑓^(𝑘)0⁺, (8) where 𝐹(𝑠)=𝑙{𝑓(𝑡);𝑠},  𝑚−1≤𝛼<𝑚.

3. System of Fractional Differential Equations and Reconstruction of Variational Iteration Method (RVIM)

Consider a system of fractional differential equations as follows: 𝐷^𝛼_∗𝑥₁(𝑡)=𝑁₁𝑡,𝑥₁,𝑥₂,…,𝑥_𝑛,𝐷^𝛼_∗𝑥₂(𝑡)=𝑁₂𝑡,𝑥₁,𝑥₂,…,𝑥_𝑛,⋮𝐷^𝛼_∗𝑥_𝑛(𝑡)=𝑁_𝑛𝑡,𝑥₁,𝑥₂,…,𝑥_𝑛, (9) where 𝑁_𝑖’s are linear/nonlinear functions of 𝑡,𝑥₁,𝑥₂,…,𝑥_𝑛,  𝐷^𝛼𝑖_∗ is the derivative of 𝑥_𝑖 with order of 𝛼_𝑖 in the sense of Caputo and 𝑚−1≤𝛼_𝑖<𝑚 with 𝑚≥1, subjected to the initial conditions: 𝑥_𝑖^(𝑘)=𝑐^𝑖_𝑘0≤𝑘≤𝛼_𝑖,1≤𝑖≤𝑛. (10) Equation (10) can be rewritten down as a correction function in the following way: 𝐷^𝛼_∗𝑥_𝑖(𝑡)=𝑁_𝑖𝑡,𝑥₁,…,𝑥_𝑛,𝑖=1,…,𝑛. (11) Therefore, the approximate solution can be reached as follows: 𝑥_𝑖(𝑡)=lim_𝑛→∞𝑥^𝑛_𝑖(𝑡),𝑖=1,…𝑛, (12) where 𝑥^𝑛_𝑖 indicates 𝑛th approximation of 𝑥_𝑖 and 𝑥⁰_𝑖 is 𝑥_𝑖(0)+𝑡𝑥^′_𝑖(0)+⋯+𝑡^𝑛𝑥_𝑖^(𝑛)(0)/𝑛! where 𝑥_𝑖^(𝑘)(0), 𝑘=0,1,…,𝑛 are substituted by initial condition of the main problem.

Taking Laplace transform to both sides of (11) in the usual way and using the homogenous initial conditions (i.e., artificial initial conditions equal to zero), the result can be obtained as follows: 𝑠^𝛼𝑖𝑙𝑥_𝑖(𝑡)=𝑙𝑁_𝑖𝑡,𝑥₁,…,𝑥_𝑛. (13) Now by applying the inverse Laplace transform to both sides of (13) and using the convolution theorem, the following relation can be concluded: 𝑥_𝑖(𝑡)=^𝑡₀(𝑡−𝜏)^𝛼−1Γ(𝛼)𝑁_𝑖𝜏,𝑥₁(𝜏),𝑥₂(𝜏),…,𝑥_𝑛(𝜏)𝑑𝜏. (14) Therefore, 𝑥_𝑖^𝑛+1(𝑡)=𝑥⁰_𝑖(𝑡)+^𝑡₀(𝑡−𝜏)^𝛼−1Γ(𝛼)×𝑁_𝑖𝜏,𝑥^𝑛₁(𝜏),𝑥^𝑛₂(𝜏),…,𝑥^𝑛_𝑛(𝜏)𝑑𝜏,𝑖=1,…,𝑛. (15) After identifying the initial approximation of 𝑥⁰_𝑖, the remaining approximations 𝑥^𝑛₁(𝑡),  𝑛>0 can be obtained. So that each term can be determined by previous term and the approximation of iteration formula can be entirely evaluated. Consequently, the solution may be written as 𝑥_𝑖(𝑡)=lim_𝑛→∞𝑥^𝑛_𝑖(𝑡),𝑖=1,…𝑛. (16)

4. Numerical Results

To demonstrate the effectiveness of the method we consider some systems of linear and nonlinear fractional differential equations.

Example 1. Let us consider the following system of two linear fractional differential equations: 𝐷^𝛽_∗𝑥(𝑡)(𝑡)(𝑡),𝐷=𝑥+𝑦_𝛾∗𝑦(𝑡)=−𝑥(𝑡)+𝑦(𝑡), (17) subjected to the initial conditions 𝑥(0)=0,𝑦(0)=1. (18) Applying the RVIM to (17), the result is as follows: 𝑙{𝑥(𝑡=1)}𝑠^𝛽𝑙𝑥(𝑡)(𝑡,𝑙𝑦(𝑡=1+𝑦))𝑠^𝛾𝑙−𝑥(𝑡)+𝑦(𝑡). (19) Applying the inverse Laplace transform to both sides of (19) results in 𝑥(𝑡)=1Γ𝛽^𝑡₀(𝑡−𝜏)^𝛽−1𝑥(𝜏)(𝜏𝑦(𝑡)=1+𝑦)𝑑𝜏,Γ𝛾^𝑡₀(𝑡−𝜏)^𝛾−1−𝑥(𝜏)+𝑦(𝜏)𝑑𝜏. (20) Therefore, approximate solution for (20) can be readily obtained as 𝑥_𝑛+1(𝑡)=𝑥₀(𝑡)+1Γ𝛽^𝑡₀(𝑡−𝜏)^𝛽−1𝑥_𝑛(𝜏)+𝑦_𝑛(𝜏𝑦)𝑑𝜏,_𝑛+1(𝑡)=𝑦₀(𝑡)+1Γ𝛾^𝑡₀(𝑡−𝜏)^𝛾−1−𝑥_𝑛(𝜏)+𝑦_𝑛(𝜏)𝑑𝜏, (21) where 𝑥₀(𝑡)=0,  𝑦₀(𝑡)=1.

According to (21), after some simplification and substitution, the following sets of equations are concluded: 𝑥₁(𝑡)=𝑡^𝛽Γ,𝑦𝛽+1₁(𝑡)𝑡=1+^𝛾Γ,𝑥𝛾+1₂(𝑡)=𝑡^𝛽Γ+𝑡𝛽+1Γ^2𝛽+𝑡2𝛽+1Γ^𝛾+𝛽,𝑦𝛾+𝛽+1₂(𝑡)𝑡=1+^𝛾Γ+𝑡𝛾+1Γ^2𝛾−𝑡2𝛾+1Γ^𝛾+𝛽,𝑥𝛾+𝛽+1₃(𝑡)=𝑡^𝛽Γ+𝑡𝛽+1Γ^2𝛽+𝑡2𝛽+1Γ^𝛾+𝛽+𝑡𝛾+𝛽+1_3𝛽Γ+𝑡3𝛽+1Γ^2𝛾+𝛽,𝑦2𝛾+𝛽+1₃(𝑡)𝑡=1+^𝛾Γ+𝑡𝛾+1Γ^2𝛾−𝑡2𝛾+1Γ^𝛾+𝛽+𝑡𝛾+𝛽+1_3𝛾Γ−𝑡3𝛾+1Γ^𝛾+2𝛽−𝛾+2𝛽+12𝑡^2𝛾+𝛽Γ.2𝛾+𝛽+1 (22) Figure 1 shows the approximate solution for system (17), obtained for the values of 𝛽=𝛾=1. This is the only case for which we know the exact solution (𝑥(𝑡)=𝑒^𝑡sin(𝑡),  𝑦(𝑡)=𝑒^𝑡cos(𝑡)). One can see that our approximate solutions by using the RVIM are in a good agreement with its exact solution.

Figure 1 

Comparison plot of RVIM method for system (17), when 𝛽=𝛾=1  with its exact solution.

Figure 2 shows the approximate solutions for system (17), obtained for the values of 𝛽=0.7  and  𝛾=0.9. It is to be noted that the following three iterations were used in evaluating the approximate solution (whereas by the differential transform method twenty-five terms were used in evaluating the approximate solutions).

(a)

(b)

Figure 2 

Plots of system (17), when 𝛽=0.7,  𝛾=0.9: (a) RVIM method and (b) differential transform method.

The results in Figures 1 and 2 are in full agreement with the results obtained in [9], using differential transform method.

Example 2. Consider the following nonlinear system of fractional differential equations: 𝑑^𝛼𝑦₁𝑑𝑡^𝛼=−𝑦₁,𝑑^𝛼𝑦₂𝑑𝑡^𝛼=𝑦₁−𝑦²₂,𝑑^𝛼𝑦₃𝑑𝑡^𝛼=𝑦²₂, (23) subjected to the initial conditions 𝑦₁(0)=1,𝑦₂(0)=0,𝑦₃(0)=0. (24) Applying the RVIM to (23), the result is as follows: 𝑙𝑦₁(𝑡=1)𝑠^𝛼𝑙−𝑦₁(𝑡,𝑙𝑦₂)(𝑡=1)𝑠^𝛼𝑙𝑦₁(𝑡)−𝑦²₂(𝑡,𝑙𝑦₃)(𝑡=1)𝑠^𝛼𝑙𝑦²₂(𝑡). (25) Benefiting from the inverse Laplace transform to both sides of (25), one obtains 𝑦₁(𝑡)=1Γ(𝛼)^𝑡₀(𝑡−𝜏)^𝛼−1−𝑦₁(𝜏𝑦)𝑑𝜏,₂(𝑡)=1Γ(𝛼)^𝑡₀(𝑡−𝜏)^𝛼−1𝑦₁(𝜏)−𝑦²₂(𝜏𝑦)𝑑𝜏,₃(𝑡)=1Γ(𝛼)^𝑡₀(𝑡−𝜏)^𝛼−1𝑦²₂(𝜏)𝑑𝜏. (26) Therefore, approximate solution for (26) can be readily obtained as 𝑦₁^𝑛+1(𝑡)=𝑦⁰₁(𝑡)+1Γ(𝛼)^𝑡₀(𝑡−𝜏)^𝛼−1−𝑦₁(𝜏𝑦)𝑑𝜏,_2𝑛+1(𝑡)=𝑦⁰₂(𝑡)+1Γ(𝛼)^𝑡₀(𝑡−𝜏)^𝛼−1𝑦₁(𝜏)−𝑦²₂(𝜏𝑦)𝑑𝜏,_3𝑛+1(𝑡)=𝑦⁰₃(𝑡)+1Γ(𝛼)^𝑡₀(𝑡−𝜏)^𝛼−1𝑦²₂(𝜏)𝑑𝜏, (27) where 𝑦⁰₁(𝑡)=1,  𝑦⁰₂(𝑡)=0,  𝑦⁰₃(𝑡)=0, and 𝑦^𝑛_𝑖 indicates 𝑛th approximation of 𝑦_𝑖 for 𝑖=1,2,3.

According to (27), after some simplification and substitution, the following set of equations is concluded: 𝑦¹₁(𝑡)𝑡=1−^𝛼Γ(),𝑦𝛼+1₁₂(𝑡)=𝑡^𝛼Γ(),𝑦𝛼+1₁₃(𝑡)𝑦=0,₂₁(𝑡)𝑡=1−^𝛼Γ()+𝑡𝛼+1Γ(^2𝛼),𝑦2𝛼+1₂₂(𝑡)=𝑡^𝛼Γ()−𝑡𝛼+1Γ(^2𝛼)−Γ()𝑡2𝛼+12𝛼+1^3𝛼(Γ(𝛼+1))²(),𝑦⋅Γ3𝛼+1₂₃(𝑡)=Γ(2𝛼+1)𝑡^3𝛼(Γ(𝛼+1))²(),𝑦⋅Γ3𝛼+1₃₁(𝑡)𝑡=1−^𝛼Γ()+𝑡𝛼+1Γ(^2𝛼),𝑦2𝛼+1₃₂(𝑡)=𝑡^𝛼Γ()−𝑡𝛼+1Γ(^2𝛼)−Γ()𝑡2𝛼+12𝛼+1^3𝛼(Γ(𝛼+1))²(),𝑦⋅Γ3𝛼+1₃₃(𝑡)=Γ(2𝛼+1)𝑡^3𝛼(Γ(𝛼+1))²().⋅Γ3𝛼+1 (28) It is to be noted that we reached the approximate solution after three iterations by the method of RVIM, whereas it is obtained after seventy iterations by the differential transform method. In Figure 3, we draw the curves of approximate solutions 𝑦₁(𝑡),  𝑦₂(𝑡), and 𝑦₃(𝑡), which is obtained for the value of 𝛼=0.7. The graphical results are in a very good agreement with the results in [9].

(a)

(b)

Figure 3 

Plots of system (23): (a) RVIM method and (b) differential transform method.

Example 3. Lastly we consider the following system of two nonlinear fractional differential equations: 𝐷_∗^1.3𝑥(𝑡)=𝑥(𝑡)+𝑦²(𝑡),𝐷_∗^2.4𝑦(𝑡)=𝑥(𝑡)+5𝑦(𝑡), (29) with the initial conditions 𝑥(0)=0,𝑥^′(0)𝑦(0)=1,=0,𝑦^′(0)=1,𝑦^′′(0)=1. (30) Applying the RVIM to (29), the result is as follows: 𝑙{𝑥(𝑡=1)}𝑠^1.3𝑙𝑥(𝑡)+𝑦²(𝑡,𝑙𝑦(𝑡=1))𝑠^2.4𝑙𝑥(𝑡)+5𝑦(𝑡). (31) Benefiting from the inverse Laplace transform to both sides of (31), one obtains 𝑥(𝑡)=1Γ(1.3)^𝑡₀(𝑡−𝜏)^0.3𝑥(𝜏)+𝑦²(𝜏𝑦(𝑡)=1)𝑑𝜏,Γ(2.4)^𝑡₀(𝑡−𝜏)^1.4𝑥(𝜏)+5𝑦(𝜏)𝑑𝜏. (32) Therefore, approximate solution for (32) can be readily obtained as: 𝑥_𝑛+1(𝑡)=𝑥₀(𝑡)+1Γ(1.3)^𝑡₀(𝑡−𝜏)^0.3𝑥(𝜏)+𝑦²(𝜏𝑦)𝑑𝜏,_𝑛+1(𝑡)=𝑦₀(𝑡)+1Γ(2.4)^𝑡₀(𝑡−𝜏)^1.4𝑥(𝜏)+5𝑦(𝜏)𝑑𝜏, (33) where the initial approximation must be satisfied by the following equations: 𝑥₀(𝑡)=𝑡,𝑦₀(𝑡)𝑡=𝑡+²2. (34) According to (33), after some simplification and substitution, the following sets of equations are concluded: 𝑥₁(𝑡)𝑡=𝑡+^2.3Γ()+3.32𝑡^3.3Γ()+4.36𝑡^5.3Γ()+6.36𝑡^4.3Γ(),𝑦5.3₁(𝑡)𝑡=𝑡+²2+6𝑡^3.4Γ()+4.45𝑡^4.4Γ().5.4 (35) Figure 4 shows the efficiency of this method to obtain approximate solutions of system (29).

(a)

(b)

Figure 4 

Plots of system (29): (a) RVIM method and (b) differential transform method.

5. Conclusion

In this paper, Reconstruction of Variational Iteration Method (RVIM) was successfully employed to solve systems of differential equations of fractional order. The work emphasized our belief that the method is a reliable technique to handle linear and nonlinear systems of fractional differential equations.

The results of this method are in a good agreement with those obtained by using the differential transform method and the Adomian decomposition method. One of the advantages of this method in comparison with the Adomian decomposition method is that we do not need to do the difficult computation for finding the Adomian polynomials.

Moreover, the method presented rapidly convergent successive approximations without any restrictive assumptions or transformation which may change the physical behavior of the problem.

Evidently, the RVIM reduced the size of calculation and also the iteration was direct and straight forward.

Generally, the proposed method is promising and applicable to a board class of linear and nonlinear systems in the theory of fractional calculus.