Mathematical Problems in Engineering
Volume 2017 (2017), Article ID 7853839, 13 pages
https://doi.org/10.1155/2017/7853839

Research Article

New Solutions for System of Fractional Integro-Differential Equations and Abel’s Integral Equations by Chebyshev Spectral Method

Hassan A. Zedan,¹ Seham Sh. Tantawy,² and Yara M. Sayed²

¹Department of Mathematics, Faculty of Science, Kafrelsheikh University, Kafrelsheikh, Egypt
²Department of Mathematics, Faculty of Education, Ain Shams University, Cairo, Egypt

Correspondence should be addressed to Hassan A. Zedan; hassanzedan2003@yahoo.com

Received 14 September 2016; Accepted 14 February 2017; Published 30 March 2017

Academic Editor: J.-C. Cortés

Copyright © 2017 Hassan A. Zedan et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Chebyshev spectral method based on operational matrix is applied to both systems of fractional integro-differential equations and Abel’s integral equations. Some test problems, for which the exact solution is known, are considered. Numerical results with comparisons are made to confirm the reliability of the method. Chebyshev spectral method may be considered as alternative and efficient technique for finding the approximation of system of fractional integro-differential equations and Abel’s integral equations.

1. Introduction

In recent years, the topic of fractional calculus has attracted many scientists because of its several applications in many areas, such as physics, chemistry, biology, and engineering. For a detailed survey with collections of applications in various fields, see, for example, [1–6].

The numerical solution of differential equations of integer order has been a hot topic in numerical and computational mathematics for a long time. There are many different methods and different basis functions have been used to estimate the solution of fractional integro-differential equations or Abel’s integral equations, such as Adomian decomposition method [7, 8], fractional differential transform method [9, 10], collocation method [11, 12], homotopy perturbation method [13, 14], homotopy analysis method [15, 16], variational iteration method [17], discrete Galerkin method [18], and Haar wavelet method [19].

Spectral methods provide a computational approach that has achieved substantial popularity over the last four decades. They have gained new popularity in automatic computations for a wide class of physical problems in fluid and heat flow. Their fascinating merit is the high accuracy. So, they have been applied successfully to numerical simulations of many problems in science and engineering; see [20–24].

The operational matrix of fractional derivatives has been determined for some types of orthogonal polynomials, such as Chebyshev polynomials [25] and Legendre polynomials [26], and for integration has been determined for several types of orthogonal polynomials, such as Chebyshev polynomials [27], Laguerre series [28], and Legendre polynomials [29]. Recently, the Bernstein operational matrix approach is developed for solving a system of high order linear Volterra–Fredholm integro-differential equations in [30].

In the present paper, we use Chebyshev spectral method based on operational matrix to solve system of fractional integro-differential equations:

with initial conditions

where .

And we use Abel’s integral equation:

where or , is a continuous function, and is constant.

2. Basic Definitions

In this section, we summarize some basic definitions and properties of fractional calculus theory.

Definition 1. A real function , is said to be in the space , if there exists a real number , such that , where Clearly if .

Definition 2. A function , is said to be in space , if .

Definition 3. The Riemann-Liouville fractional integral operator of order , of a function , is defined as

Definition 4. The Caputo fractional derivatives of order are defined as

where and is the classical differential operator of order

For Caputo derivative, we have

We use the ceiling function denoting the smallest integer greater than or equal to and the floor function denoting the largest integer less than or equal to Also and . Recall that, for , the Caputo differential operator coincides with the usual differential operator of integer order.

More properties of the fractional derivatives and the fractional integral can be found in [3, 4].

3. Some Properties of the Shifted Chebyshev Polynomials

The well-known Chebyshev polynomials are defined on the interval and can be determined with the aid of the following recurrence formula:

where and The Chebyshev polynomials are orthogonal on the interval with respect to the weight function These polynomials satisfy the relation

where

The analytic form of the Chebyshev polynomial of degree is given by

where denotes the integer part of The zeros of are denoted by

In order to use these polynomials on the interval , we defined the so-called shifted Chebyshev polynomials by introducing the change of variable Let the shifted Chebyshev polynomials be denoted by , satisfying the orthogonality relation

where .

The shifted Chebyshev polynomials are defined as and the analytic form is given by

where and

In this form, may be generated with the aid of the following recurrence formula:

where and The zeros of are denoted by

A function , square integrable in , may be expressed in terms of the shifted Chebyshev polynomials as

where the coefficients are given by

In practice, only the first -terms shifted Chebyshev polynomials are considered. Hence, if we write

where the shifted Chebyshev coefficient vector and the shifted Chebyshev vector are given by

then the derivative of the vector can be expressed by

where is the operational matrix of derivative given by

for example, for even , we have

4. The Shifted Chebyshev Operational Matrix (COM) Fractional Derivatives

The main objective of this subsection is to generalize the COM of derivatives for the fractional calculus. By using (19), it is clear that

where and the superscript, in , denotes matrix powers. Thus

Lemma 5. Let be a shifted Chebyshev polynomial; then

Theorem 6. Let be the shifted Chebyshev vector defined in (18) and suppose ; then

where is the COM of derivatives of order in the Caputo sense and is defined as follows:

where

Note that, in , the first rows are all zero.

Proof. See [25].

Remark 7. If , then Theorem 6 gives the same result as (22).

5. System of Fractional Integro-Differential Equations

In order to use COM for system of fractional integro-differential equations of the form (1), we first approximate , and by the shifted Chebyshev polynomials as

By substituting these equations in (1), we get

Also, by substituting (22) and (28) in (2), we obtain

Then we have to collocate (31) at the shifted Chebyshev roots . These equations together with (32) generate of algebraic equations which can then be solved for the unknown coefficients of the vectors , using a suitable method. Consequently, the approximate solution , can be obtained.

In our computations we used the Gaussian elimination method to solve the resulting linear system of algebraic equations and Newton’s iteration method to solve the resulting nonlinear system of algebraic equations.

Now, we can present the following problems.

Example 8. Consider the following system of fractional integro-differential equations [7, 15]:

with the initial conditions

The exact solutions, when , are

We use Chebyshev spectral method; we may write the approximate solution

where

Substituting (36) in (33) and (34), for , we get

The roots of the shifted Chebyshev polynomial are given by

Now, for calculating the shifted Chebyshev coefficient for , substitute (40) in (38), and solving the resulting linear system of equations and (39), we get

Thus using (36), we get

Table 1 shows the comparison between the exact solution and the approximate solution with the absolute error at , and Table 2 shows the maximum of absolute error between exact solutions and approximate solutions for various choices of Figure 1 shows the graph of the exact solution and the approximate solution at , . Figure 2 shows the graph of the exact solutions and the approximate solutions at , , and .

Table 1

1
1
0
−1
−1
6.66134 × 10⁻¹⁶

1.20517
1.20517
5.23892 × 10⁻¹¹
−1.00517
−1.00517
4.22691 × 10⁻¹¹

1.42140
1.42140
4.84666 × 10⁻¹¹
−1.02140
−1.02140
5.00378 × 10⁻¹¹

1.64986
1.64986
1.00231 × 10⁻¹²
−1.04986
−1.04986
4.82325 × 10⁻¹²

1.89182
1.89182
9.13603 × 10⁻¹¹
−1.09182
−1.09182
8.16736 × 10⁻¹¹

2.14872
2.14872
1.82374 × 10⁻¹¹
−1.14872
−1.14872
1.17573 × 10⁻¹²

2.42212
2.42212
6.51625 × 10⁻¹¹
−1.22212
−1.22212
7.4986 × 10⁻¹¹

2.71375
2.71375
1.55334 × 10⁻¹¹
−1.31375
−1.31375
1.77258 × 10⁻¹¹

3.02554
3.02554
7.43441 × 10⁻¹¹
−1.42554
−1.42554
6.87448 × 10⁻¹¹

3.35960
3.35960
2.49547 × 10⁻¹¹
−1.55960
−1.55960
4.03237 × 10⁻¹¹

3.71828
3.71828
2.28355 × 10⁻¹¹
−1.71828
−1.71828
1.71747 × 10⁻¹¹

Table 2

Maximum of absolute error of
Maximum of absolute error of

2
5.35939 × 10⁻²
3.73645 × 10⁻²

4
9.25152 × 10⁻⁵
6.71757 × 10⁻⁵

6
7.59171 × 10⁻⁷
5.37812 × 10⁻⁷

8
9.13603 × 10⁻¹¹
8.16736 × 10⁻¹¹

(a)

(b)

Figure 1: The graphs of the exact solutions and the approximate solutions at , .

(a)

(b)

Figure 2: The graphs of the exact solutions and the approximate solutions at , , and 0.75.

Example 9. Consider the following nonlinear fractional system of integro-differential equations [15]:

with the initial conditions

The exact solutions, when , are

We use Chebyshev spectral method; we may write the approximate solution

where

Substituting (46) in (43) and (44), for , we get

The roots of the shifted Chebyshev polynomial are

Now, for calculating the shifted Chebyshev coefficient for , substitute (50) in (48), and solving the resulting nonlinear system of equations and (49), we get

Thus using (46), we get

Table 3 shows the comparison between the exact solution and the approximate solution with the absolute error at . Figure 3 shows the graph of the exact solution and the approximate solution at , . Figure 4 shows the graph of the exact solutions and the approximate solutions at , , and .

Table 3

1
2.22045 × 10⁻¹⁶
−1
−1

1.2048868
2.84157 × 10⁻⁴
−1.0051709
−1.0048309
3.39987 × 10⁻⁴

1.4206154
7.87348 × 10⁻⁴
−1.0214028
−1.0204521
9.50662 × 10⁻⁴

1.6487885
1.07031 × 10⁻³
−1.0498588
−1.0485561
1.30275 × 10⁻³

1.8910086
8.16141 × 10⁻⁴
−1.0918247
−1.0908354
9.89323 × 10⁻⁴

2.1488781
1.56876 × 10⁻⁴
−1.1487213
−1.1489826
2.61337 × 10⁻⁴

2.4239998
1.88101 × 10⁻³
−1.2221188
−1.2246903
2.57152 × 10⁻³

2.7179761
4.22339 × 10⁻³
−1.3137527
−1.3196511
5.89836 × 10⁻³

3.0324095
6.86862 × 10⁻³
−1.42554092
−1.4355574
1.00165 × 10⁻²

3.3689027
9.2996 × 10⁻³
−1.5596031
−1.5741019
1.44988 × 10⁻²

3.7290581
1.07763 × 10⁻²
−1.7182818
−1.7369771
1.86953 × 10⁻²

(a)

(b)

Figure 3: The graphs of the exact solutions and the approximate solutions at , .

(a)

(b)

Figure 4: The graphs of the exact solutions and the approximate solutions at , , 1.8, and 1.7.

6. Abel’s Integral Equation

In order to use COM for Abel’s integral equation of the form (3), we first approximate by the shifted Chebyshev polynomials as

By substituting (53) in (3), we get

Then we have to collocate (54) at the shifted Chebyshev roots . These equations generate linear algebraic equations which can be solved for the unknown coefficients of the vector , using a suitable method. Consequently, given in (53) can be calculated, which gives a solution of (3).

Example 10. Consider Abel’s integral equation of the first kind [31, 32]

which has the exact solution .

By applying the Chebyshev spectral method, we may write the approximate solution

Substituting (56) in (55), we get

The roots of the shifted Chebyshev polynomial are

Now, calculating the shifted Chebyshev coefficient for by substituting (58) in (57) and solving four equations yields

Therefore, we have

which is the exact solution.

Example 11. Consider Abel’s integral equation of the first kind [32]

which has the exact solution

Table 4 shows the comparison between the exact solution and the approximate solution with the absolute error at , and Table 5 shows the maximum of absolute error between exact solution and approximate solution for various choices of Figure 5 shows the graph of the exact solution and the approximate solution at .

Table 4

0
0
−1.9732143095 × 10⁻⁵
1.97321 × 10⁻⁵

0.1
0.0316227766
0.03162341809
6.41487 × 10⁻⁷

0.2
0.0894427191
0.08944237293
3.46169 × 10⁻⁷

0.3
0.1643167673
0.16431491828
1.84897 × 10⁻⁶

0.4
0.2529822128
0.25298265879
4.45981 × 10⁻⁷

0.5
0.3535533906
0.35355438911
9.98519 × 10⁻⁷

0.6
0.4647580015
0.46475676796
1.23323 × 10⁻⁶

0.7
0.5856620186
0.58566199942
1.77613 × 10⁻⁸

0.8
0.7155417528
0.71554297116
1.21836 × 10⁻⁶

0.9
0.8538149682
0.85381621495
1.2467 × 10⁻⁶

1
1
1.00000610948
6.10948 × 10⁻⁶

Table 5

Maximum of absolute error

4
1.90301 × 10⁻³

8
2.57352 × 10⁻⁴

12
8.39682 × 10⁻⁵

16
3.73104 × 10⁻⁵

20
1.97321 × 10⁻⁵

Figure 5: The graph of the exact solution and the approximate solution at .

Example 12. Consider Abel’s integral equation of the second kind [32, 33]

which has the exact solution

By applying the Chebyshev spectral method, we may write the approximate solution

Substituting (63) in (62), we get

The roots of the shifted Chebyshev polynomial are

Now, calculating the shifted Chebyshev coefficient for by substituting (65) in (64) and solving three equations yields

Therefore, we have

which is the exact solution.

Example 13. Consider Abel’s integral equation of the second kind [32, 34]

which has the exact solution

We use Chebyshev spectral method; we may write the approximate solution

Substituting (69) in (68), we get

The roots of the shifted Chebyshev polynomial are

Now, calculating the shifted Chebyshev coefficient for by substituting (71) in (70) and solving three equations yields

Therefore, we have

which is the exact solution.

Example 14. Consider Abel’s integral equation of the second kind [32, 34]

which has the exact solution

Table 6 shows the comparison between the exact solution and the approximate solution with the absolute error at , and Table 7 shows the maximum of absolute error between exact solution and approximate solution for various choices of Figure 6 shows the graph of the exact solution and the approximate solution at .

Table 6

0
1
1
2.22045 × 10⁻¹⁶

0.1
0.9090909091
0.9090909091
2.22045 × 10⁻¹⁶

0.2
0.8333333333
0.8333333333
1.77636 × 10⁻¹⁵

0.3
0.7692307692
0.76923076923
6.66134 × 10⁻¹⁶

0.4
0.7142857143
0.7142857143
3.33067 × 10⁻¹⁶

0.5
0.6666666667
0.6666666667
3.55271 × 10⁻¹⁵

0.6
0.625
0.625
2.22045 × 10⁻¹⁶

0.7
0.588235294
0.588235294
6.66134 × 10⁻¹⁶

0.8
0.555555556
0.555555556
1.11022 × 10⁻¹⁶

0.9
0.526315790
0.526315790
1.55431 × 10⁻¹⁵

1
0.5
0.5
2.44249 × 10⁻¹⁵

Table 7

Maximum of absolute error

4
2.47780 × 10⁻⁴

8
2.05220 × 10⁻⁷

12
1.80403 × 10⁻¹⁰

16
1.55986 × 10⁻¹³

20
3.55271 × 10⁻¹⁵

Figure 6: The graph of the exact solution and the approximate solution at .

Example 15. Consider the linear system of singular Volterra integral equations [10]

The exact solutions are

Table 8 shows the comparison between the exact solution and the approximate solution with the absolute error at , and Table 9 shows the maximum of absolute error between exact solution and approximate solution for various choices of Figure 7 shows the graph of the exact solution and the approximate solution at .

Table 8

1.62145 × 10⁻⁴
1.62145 × 10⁻⁴
0
2.38416 × 10⁻²
2.38416 × 10⁻²

0.0999982
1.76511 × 10⁻⁶
0.316228
0.316404
1.7626 × 10⁻⁴

0.19999
9.70125 × 10⁻⁶
0.447214
0.447335
1.21557 × 10⁻⁴

0.299981
1.92749 × 10⁻⁵
0.547723
0.547788
6.49502 × 10⁻⁵

0.399991
9.32098 × 10⁻⁶
0.632456
0.632369
8.61806 × 10⁻⁵

0.49999
9.57004 × 10⁻⁶
0.707107
0.70709
1.69007 × 10⁻⁵

0.599978
2.1643 × 10⁻⁵
0.774597
0.774625
2.87108 × 10⁻⁵

0.699983
1.74306 × 10⁻⁵
0.83666
0.836603
5.70565 × 10⁻⁵

0.799983
1.74232 × 10⁻⁵
0.894427
0.894364
6.33597 × 10⁻⁵

0.899978
2.16216 × 10⁻⁵
0.948683
0.948623
6.03199 × 10⁻⁵

0.999976
2.44081 × 10⁻⁵
1
0.999923
7.72476 × 10⁻⁵

Table 9

Maximum of absolute error of
Maximum of absolute error of

4
3.08564 × 10⁻³
1.02329 × 10⁻¹

8
9.08442 × 10⁻⁴
5.59578 × 10⁻²

12
4.28561 × 10⁻⁴
3.85957 × 10⁻²

16
2.486 × 10⁻⁴
2.9472 × 10⁻²

20
1.62145 × 10⁻⁴
2.38416 × 10⁻²

(a)

(b)

Figure 7: The graphs of the exact solutions and the approximate solutions at .

7. Conclusion

In this article, we develop the Chebyshev spectral method based on operational matrix for solving linear and nonlinear system of fractional integro-differential equations and Abel’s integral equations. Our approach was based on the shifted Chebyshev collocation methods. It can be concluded that the Chebyshev spectral method is very powerful and efficient technique for finding exact solutions for wide classes of problems.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

References

K. S. Miller and B. Ross, An Introduction to the Fractional Calculus and Fractional Differential Equations, John Wiley & Sons, New York, NY, USA, 1993.

R. Hilfer, Applications of Fractional Calculus in Physics, World Scientific, Singapore, 2000.

S. G. Samko, A. A. Kilbas, and O. I. Marichev, Fractional Integrals and Derivatives: Theory and Applications, Gordon and Breach, London, UK, 1993.

I. Podlubny, Fractional Differential Equations, vol. 198 of Mathematics in Science and Engineering, Academic Press, New York, NY, USA, 1999.

A. Debbouche and J. J. Nieto, “Sobolev type fractional abstract evolution equations with nonlocal conditions and optimal multi-controls,” Applied Mathematics and Computation, vol. 245, pp. 74–85, 2014.

A. Debbouche and J. J. Nieto, “Relaxation in controlled systems described by fractional integro-differential equations with nonlocal control conditions,” Electronic Journal of Differential Equations, vol. 2015, no. 89, pp. 1–18, 2015.

S. Momani and R. Qaralleh, “An efficient method for solving systems of fractional integro-differential equations,” Computers and Mathematics with Applications, vol. 52, no. 3-4, pp. 459–470, 2006.

R. H. Khan and H. O. Bakodah, “Adomian decomposition method and its modification for nonlinear Abel's integral equation,” International Journal of Mathematical Analysis, vol. 7, no. 45-48, pp. 2349–2358, 2013.

A. Arikoglu and I. Ozkol, “Solution of fractional integro-differential equations by using fractional differential transform method,” Chaos, Solitons & Fractals, vol. 40, no. 2, pp. 521–529, 2009.

H. Taghvafard and G. H. Erjaee, “On solving a system of singular Volterra integral equations of convolution type,” Communications in Nonlinear Science and Numerical Simulation, vol. 16, no. 9, pp. 3486–3492, 2011.

J. Zhao, J. Xiao, and N. J. Ford, “Collocation methods for fractional integro-differential equations with weakly singular kernels,” Numerical Algorithms, vol. 65, no. 4, pp. 723–743, 2014.

R. K. Pandey, S. Sharma, and K. Kumar, “Collocation method for generalized Abel's integral equations,” Journal of Computational and Applied Mathematics, vol. 302, pp. 118–128, 2016.

K. Sayevand, M. Fardi, E. Moradi, and F. Hemati Boroujeni, “Convergence analysis of homotopy perturbation method for Volterra integro-differential equations of fractional order,” Alexandria Engineering Journal, vol. 52, no. 4, pp. 807–812, 2013.

S. Kumar, O. P. Singh, and S. Dixit, “Homotopy perturbation method for solving system of generalized Abel's integral equations,” Applications and Applied Mathematics, vol. 6, no. 11, pp. 2009–2024, 2011.

M. Zurigat, S. Momani, and A. Alawneh, “Homotopy analysis method for systems of fractional integro-differential equations,” Neural, Parallel and Scientific Computations, vol. 17, no. 2, pp. 169–186, 2009.

F. Mirzaee, M. Komak Yari, and M. Paripour, “Solving linear and nonlinear Abel fuzzy integral equations by homotopy analysis method,” Journal of Taibah University for Science, vol. 9, no. 1, pp. 104–115, 2015.

H. A. Zedan, S. S. Tantawy, and Y. M. Sayed, “Convergence of the variational iteration method for initial-boundary value problem of fractional integro-differential equations,” Journal of Fractional Calculus and Applications, vol. 5, supplement 3, pp. 1–14, 2014.

P. Mokhtary, “Discrete Galerkin method for fractional integro-differential equations,” Acta Mathematica Scientia. Series B. English Edition, vol. 36, no. 2, pp. 560–578, 2016.

H. Saeedi, N. Mollahasani, M. Mohseni Moghadam, and G. N. Chuev, “An operational Haar wavelet method for solving fractional Volterra integral equations,” International Journal of Applied Mathematics and Computer Science, vol. 21, no. 3, pp. 535–547, 2011.

C. Canuto, M. Y. Hussaini, A. Quarteroni, and T. A. Zang, Spectral Methods in Fluid Dynamics, Springer, New York, NY, USA, 1988.

E. A. Coutsias, T. Hagstrom, and D. Torres, “An efficient spectral method for ordinary differential equations with rational function coefficients,” Mathematics of Computation, vol. 65, no. 214, pp. 611–635, 1996.

Y. Yang, Y. Chen, and Y. Huang, “Spectral-collocation method for fractional Fredholm integro-differential equations,” Journal of the Korean Mathematical Society, vol. 51, no. 1, pp. 203–224, 2014.

E. H. Doha, A. H. Bhrawy, and M. A. Saker, “Integrals of Bernstein polynomials: an application for the solution of high even-order differential equations,” Applied Mathematics Letters, vol. 24, no. 4, pp. 559–565, 2011.

E. H. Doha and A. H. Bhrawy, “Efficient spectral-Galerkin algorithms for direct solution of fourth-order differential equations using Jacobi polynomials,” Applied Numerical Mathematics, vol. 58, no. 8, pp. 1224–1244, 2008.

E. H. Doha, A. H. Bhrawy, and S. S. Ezz-Eldien, “A Chebyshev spectral method based on operational matrix for initial and boundary value problems of fractional order,” Computers & Mathematics with Applications, vol. 62, no. 5, pp. 2364–2373, 2011.

A. Saadatmandi and M. Dehghan, “A new operational matrix for solving fractional-order differential equations,” Computers and Mathematics with Applications. An International Journal, vol. 59, no. 3, pp. 1326–1336, 2010.

A. H. Bhrawy and A. S. Alofi, “The operational matrix of fractional integration for shifted Chebyshev polynomials,” Applied Mathematics Letters, vol. 26, no. 1, pp. 25–31, 2013.

C. Hwang and Y. P. Shih, “Parameter identification via Laguerre polynomials,” International Journal of Systems Science, vol. 13, no. 2, pp. 209–217, 1982.

P. N. Paraskevopoulos, “Legendre series approach to identification and analysis of linear systems,” IEEE Transactions on Automatic Control, vol. 30, no. 6, pp. 585–589, 1985.

K. Maleknejad, B. Basirat, and E. Hashemizadeh, “A Bernstein operational matrix approach for solving a system of high order linear Volterra–Fredholm integro-differential equations,” Mathematical and Computer Modelling, vol. 55, no. 3-4, pp. 1363–1372, 2012.

S. Sohrabi, “Comparison Chebyshev wavelets method with BPFs method for solving Abel's integral equation,” Ain Shams Engineering Journal, vol. 2, no. 3-4, pp. 249–254, 2011.

Z. Avazzadeh, B. Shafiee, and G. B. Loghmani, “Solution of Abels integral equations using legendre polynomials and fractional calculus techniques,” International Journal of Mathematical Archive (IJMA), vol. 2, no. 8, pp. 1352–1359, 2011.

E. A. Galperin and E. J. Kansa, “Application of global optimization and radial basis functions to numerical solutions of weakly singular Volterra integral equations,” Computers & Mathematics with Applications, vol. 43, no. 3-5, pp. 491–499, 2002.

R. K. Pandey, O. P. Singh, and V. K. Singh, “Efficient algorithms to solve singular integral equations of Abel type,” Computers & Mathematics with Applications, vol. 57, no. 4, pp. 664–676, 2009.