Subscribe Now Subscribe Today
Research Article

Numerical Solution of Fractional Ordinary Differential Equations using a Multiquadric Approximation Scheme

M. Avaji, J.S. Hafshejani, S.S. Dehcheshmeh and D.F. Ghahfarokhi
Facebook Twitter Digg Reddit Linkedin StumbleUpon E-mail

In this study, the aim was to present a Multiquadric (MQ) approximation scheme on the numerical solution of Fractional Ordinary Differential Equations (FODEs) using extensively in engineering. The properties of the proposed MQ approximation scheme and its advantages which include using data points in arbitrary locations with arbitrary ordering are presented. Comparing between the numerical results obtained from our method and the other methods confirms the good accuracy of the presented scheme. Also, we present some experiments wherein the numerical results show that the method works excellently, even where the data points are scattered.

Related Articles in ASCI
Similar Articles in this Journal
Search in Google Scholar
View Citation
Report Citation

  How to cite this article:

M. Avaji, J.S. Hafshejani, S.S. Dehcheshmeh and D.F. Ghahfarokhi, 2012. Numerical Solution of Fractional Ordinary Differential Equations using a Multiquadric Approximation Scheme. Journal of Applied Sciences, 12: 168-173.

DOI: 10.3923/jas.2012.168.173

Received: August 03, 2011; Accepted: November 22, 2011; Published: January 18, 2012


The theory of fractional calculus was first raised in the year 1695 by Marquis de L’ Hopital and from now on many studies were done and many important books were published in this field wherein we can point out to the study of Oldham and Spanier (1974), Miller and Ross (1993), Samko et al. (1993) and Podlubny (1999). Most of the scientific problems and phenomena are modeled by Fractional Ordinary Differential Equations (FODEs) and Fractional Partial Differential Equations (FPDEs). For instance; in mathematical physics (Podlubny, 1999), in fluid and continuum mechanics (Carpinteri and Mainardi, 1997), coloured noises (Sun et al., 1984), biology, chemistry, acoustics and psychology (Ahmad and El-Khazali, 2007). Some of FPDEs have been studied and solved including the classic Fractional PDEs (Vanani and Aminataei, 2011a), the fractional KdV equation (Momani, 2005) and linear and nonlinear space and time-fractional diffusion-wave equation (Momani et al., 2007; Jafari and Seifi, 2009).

In most cases, these problems do not admit analytical solution, so these equations should be solved using special techniques. In the last decade, several computational methods have been applied to solve FDEs, prominent among which are the Homotopy Perturbation Method (HPM) (Momani and Odibat, 2007a; Jafari and Seifi, 2009), the Adomian Decomposition Method (ADM) (El-Sayed and Gaber, 2006; Momani and Odibat, 2006), the Variational Iteration Method (VIM) (Momani and Odibat, 2007b), the Generalized Differential Transformation Method (GDTM) (Momani and Odibat, 2008) and the Fractional Difference Method (FDM) (Momani and Odibat, 2007a; Ghorbani, 2008).

Out of these methods, we would like to present a simple and efficient method for solving FODEs. The MQ approximation scheme is an useful method for the numerical solution of ordinary and partial differential equations (ODEs and PDEs). It is a grid-free spatial approximation scheme that converges exponentially for the spatial terms of ODEs and PDEs.

The MQ approximation scheme was first introduced by Hardy (1971) who successfully applied this method to approximate surfaces and bodies from field data. Hardy (1990) has written a detailed review article summarizing its explosive growth since it was first introduced. In 1972, Franke (1982) published a detailed comparison of 29 different scattered data schemes for analytic problems. Of all the techniques tested, he concluded that MQ performed the best in accuracy, visual appeal and ease of implementation, even against various finite element schemes. To the best of our knowledge, this is the first demonstration of the application of the MQ approximation scheme to FODEs.


Here, we state some preliminaries and definitions of fractional calculus (Podlubny, 1999):

Definition 1: A real function u (x), x>0 is said to be in the space Cμ, μεR, if there exists a real number p>μ such that u (x) = xpv (x), where v(x)εC (0, ∞) and it is said to be in the space Cμm iff u(m) (x)εCμ, mεN.

Definition 2: The Riemann--Liouville fractional integral operator of order α≥0, of a function u (x)εCμ, μ≥-1, is defined as:

where, Γ is the Gamma function. Some of the most important properties of operator Jα for u (x)εCμ, μ≥-1, α, β≥0 and γ>-1, are as follows:

Jα Jβu (x) = J(α+β) u (x)
Jα Jβu (x) = Jβ Jα u (x)

Definition 3: The fractional derivative of u (x) in the Caputo’s sense is defined as:


where, m-1<α≤m, mεN, x>0, u (x)εCm-1.


The basic MQ approximation scheme assumes that any function can be expanded as a finite series of upper hyperboloids that are written as follows:


where, N is the total number of data centers under consideration and:


(x-xj)2 is the square of Euclidean distances in R and R2>0 is an input shape parameter. Note that, the basis function h, which is a type of spline approximation, is continuously differentiable. The expansion coefficients aj are found by solving a set of full linear equations, which are written as follows:


Zerroukut et al. (1998) found that a constant shape parameter (R2) achieves better accuracy. Mai-Duy and Tran-Cong (2001) have developed new methods based on Radial Basis Function Networks (RBFNs) to approximate both functions and their first and higher derivatives. The Direct RBFN (DRBFN) and Indirect RBFN (IRBFN) methods have been studied and it has been found that the IRBFN methods yield consistently better results for both functions and their derivatives. Recently, Aminataei and Mazarei (2005) stated that, in the numerical solution of elliptic PDEs using DRBFN and IRBFN methods, the IRBFN method is more accurate than other methods with very small error. They have shown that, especially, on one-dimensional equations, IRBFN method is more accurate than DRBFN method. Furthermore, Aminataei and Mazarei (2008) used the DRBFN and IRBFN methods on the polar coordinate and have achieved better accuracy.

Micchelli (1986) proved that MQ schemes belong to a class of conditionally positive definite RBFNs. He showed that equation is always solvable for distinct points. Madych and Nelson (1990) proved that MQ interpolation always produces a minimal semi-norm error and that the MQ interpolant and derivative estimates converge exponentially as the density of data centers increases.

In contrast, the MQ interpolant is continuously differentiable over the entire domain of data centers and the spatial derivative approximations were found to be excellent, especially in very steep gradient regions where traditional methods fail. This ability to approximate spatial derivatives is due in large part to a slight modification of the original MQ scheme that permits the shape parameter to vary with the basis function.

Instead of using the constant shape parameter defined in Eq. 3, we have used variable shape parameters (Kansa, 1990a, b; Vanani and Aminataei, 2009, 2008, 2011b) as follows:





R2max and R2min are two input parameters chosen so that the following ratio is in the given range:


Madych (1992) has proven that very large values of a shape parameter are desirable in certain circumstances. Equation 5 is one way to have at least one very large value of a shape parameter without incurring the onset of severe ill-conditioned problems.

Now, in order to apply the MQ approximation scheme for solving FODEs, let us consider a FODE in the form:

Dαu (x)+Lu (x) = f (x), 0≤x≤b

u(k) (0) = λk, k = 0,1,…, [α]

where, L is the differential operator, Dα is the fractional differential operator of order α that operates on the interior and f: [0, b]→R is a known function.

Let be the N-[α] collocation points in [0, b]. Substituting the collocation points into Eq. 8, we obtain:


Therefore, we have the following system:




Therefore, imposing the [α] initial conditions, the system of N equations with N unknowns is available. Then, we must solve this system to make distinct the unknown coefficients. Hence, we have used the Gauss elimination method with total pivoting to solve such a system.

Remark: It is noticeable that collocating points can be scattered. This is one of the most important advantages of the MQ approximation scheme [0, 0]. The numerical results show this issue easily and the applicability of the MQ approximation scheme in this sense, is observable.


Here, three experiments including linear and nonlinear problems on regular and irregular domains are solved using MQ approximation scheme. Since, most of FDEs have not exact solution, therefore we must compare them with known numerical methods such as HPM, VIM, FDM and ADM. The results and their comparison with several powerful methods illustrate the validity and capability of MQ approximation scheme. If the exact solution of the problem exists, the accuracy of an approximate solution is measured by means of the discrete relative L2 norm defined as:


where, u and are the exact and computed solutions, respectively and N is the number of unknown nodal values of u. The computations associated with the experiments were also performed in Maple 14 on a PC, CPU 2.8 GHz.

Experiment 1: Consider the following FODE:

Dα u (x)+u″ (x)+u (x) = 8, x>0, 0<α≤2

with the initial conditions: u (0) = u’ (0) = 0

This problem is chosen from Momani and Odibat (2007a) and Arikoglu and Ozkol (2007). It was solved by FDM, ADM and VIM methods. Here, we have solved it by MQ approximation scheme with Rmax = 250, Rmin = 0.5, N = 30 for different α and have compared it with the closed form series solution of the exact solution (Momani and Odibat, 2007b; Arikoglu and Ozkol, 2007) and the three aforesaid mentioned methods. The results show a good agreement with the same results obtained by using the MQ approximation scheme. The comparison is shown in Table 1.

In the above experiment, we have compared the MQ approximation scheme with the FDM, ADM and VIM for solving a FODE. We observe that higher-order accuracy can be achieved by MQ-RBF approximation scheme than using the same terms in the other methods. The discrete relative L2 norm (Ne) given by Eq. 12 for α = 0.5, 1 and 1.5 is 0.000009, 0.0 and 0.000467, respectively. These results, demonstrate another capability of the proposed method. Furthermore, when α is an integer, the MQ approximation scheme is so much accurate than the other methods.

Experiment 2: Consider the following FODE (Abdulaziz et al., 2008):

Dα u (x) = -u (x), 0<α<1, x>0, u (0) = 1

The close form series solution of the exact solution is obtained by HPM [0] as:

We have tested it for seven scattered data points (N = 7) for α = 0.75, 0.85, 0.95 with Rmax = 500, Rmin = 0.3 and have compared it with the HPM. The comparison is given in Table 2.

Table 1: Comparison of the solutions of FDM, ADM, VIM and MQ-RBF approximationscheme with the exact solution for different α and x of the experiment 1

The computed results in Table 2, demonstrate the approximate solution obtained using the MQ approximation scheme, is in good agreement with the approximate solution obtained using the HPM for all values of scattered data points.

Experiment 3: Consider the following nonlinear FODE:

Dα u (x)+eu (x) = 0, 0≤x≤1, 0<α≤1

u (0) = 0

The exact solution of this initial value problem for α = 1, is u (x) = -In (x+1).

In this experiment, we use from the following expansion to linearize eu (x) as:

We have solved this experiment for N = 11 and α = 0.75, 0.85, 0.95, 0.99, 0.999.

Figure 1, shows the numerical results and illustrates a good agreement with the same results obtained by using the VIM and ADM in [0].

Fig. 1: The graphs of the experiment 4.3 for N = 11 and different α

Table 2: Comparison of the solutions of HPM and MQ-RBF approximation scheme for different α and scattered data points of X of the experiment 2

Also, the results show that when α→1, then approximate solution tends to exact solution, rapidly.


In this study, a MQ approximation scheme is proposed to solve FODEs. The results reveal that the technique introduced here is effective and convenient in solving FODEs. This method is also easy to implement and yields the desired accuracy with only a few terms. Other advantages of the present method are its low run time, a minimal number of data points in the required domain and its applicability to scattered collocated points. All of these advantages of the MQ approximation scheme suggest that the method is a valid and powerful tool.

1:  Oldham, K.B. and J. Spanier, 1974. The Fractional Calculus: Theory and Applications of Differentiation and Integration to Arbitrary Order. Elsevier Academic Press, New York and London, ISBN-13: 9780125255509, Pages: 234.

2:  Miller, K.S. and B. Ross, 1993. An Introduction to The Fractional Calculus and Fractional Differential Equations. 1st Edn., John Wiley and Sons Inc., New York, ISBN: 0471588849, Pages: 384.

3:  Samko, S.G., A.A. Kilbas and O.I. Marichev, 1993. Fractional Integrals and Derivatives: Theory and Applications. Gordon and Breach Sciench Publishers, Switzerland.

4:  Podlubny, I., 1999. Fractional Differential Equations. 1st Edn., Academic Press, New York.

5:  Carpinteri, A. and F. Mainardi, 1997. Fractals and Fractional Calculus Continuum Mechanics. Springer, New York, USA., ISBN-13: 9783211829134, pp: 291-348.

6:  Sun, H.H., A.A. Abdelwahab and B. Onaral, 1984. Linear approximation of transfer function with a pole of fractional power. IEEE Trans. Automat. Control, 29: 441-444.
Direct Link  |  

7:  Ahmad, W.M. and R. El-Khazali, 2007. Fractional-order dynamical models of love. Chaos Solitons Fractals, 33: 1367-1375.
CrossRef  |  Direct Link  |  

8:  Vanani, S.K. and A. Aminataei, 2011. Tau approximate solution of fractional partial differential equations. Comput. Math. Appl., 62: 1075-1083.
CrossRef  |  Direct Link  |  

9:  Momani, S., 2005. An explicit and numerical solutions of the fractional KdV equation. Math. Comput. Simul., 70: 110-118.
CrossRef  |  Direct Link  |  

10:  Momani, S., Z. Odibat and V.S. Erturk, 2007. Generalized differential transform method for solving a space- and time-fractional diffusion-wave equation. Phys. Lett. A, 370: 379-387.
CrossRef  |  Direct Link  |  

11:  Jafari, H. and S. Seifi, 2009. Homotopy analysis method for solving linear and nonlinear fractional diffusion-wave equation. Commun. Nonlinear Sci. Numer. Simul., 14: 2006-2012.
CrossRef  |  Direct Link  |  

12:  Momani, S. and Z. Odibat, 2007. Comparison between the homotopy perturbation method and the variational iteration method for linear fractional partial differential equations. Comput. Math. Appl., 54: 910-919.
Direct Link  |  

13:  El-Sayed, A.M.A. and M. Gaber, 2006. The Adomian decomposition method for solving partial differential equations of fractal order in finite domains. Phys. Lett. A, 359: 175-182.
CrossRef  |  Direct Link  |  

14:  Momani, S. and Z. Odibat, 2006. Analytical solution of a time-fractional Navier-Stokes equation by Adomian decomposition method. Applied Math. Comput., 177: 488-494.
CrossRef  |  Direct Link  |  

15:  Momani, S. and Z. Odibat, 2008. A novel method for nonlinear fractional partial differential equations: Combination of DTM and generalized Taylor's formula. J. Comput. Applied Math., 220: 85-95.
CrossRef  |  Direct Link  |  

16:  Momani, S. and Z. Odibat, 2007. Numerical comparison of methods for solving linear differential equations of fractional order. Chaos Solitons Fractals, 31: 1248-1255.
CrossRef  |  

17:  Ghorbani, A., 2008. Toward a new analytical method for solving nonlinear fractional differential equations. Comput. Methods Applied Mech. Eng., 197: 4173-4179.
CrossRef  |  

18:  Hardy, R.L., 1971. Multiquadric equations of topography and other irregular surfaces. J. Geophys. Res., 76: 1905-1915.
CrossRef  |  Direct Link  |  

19:  Hardy, R.L., 1990. Theory and applications of the multiquadric bi-harmonic method: 20 years of discovery. Comput. Math. Appl., 19: 163-208.

20:  Franke, R., 1982. Scattered data interpolation: Tests of some methods. Math. Comput., 38: 181-200.
CrossRef  |  Direct Link  |  

21:  Zerroukut, M., H. Power and C.S. Chen, 1998. A numerical method for heat transfer problems using collocation and radial basis functions. Int. J. Numer. Math. Eng., 42: 1263-1278.
CrossRef  |  Direct Link  |  

22:  Mai-Duy, N. and T. Tran-Cong, 2001. Numerical solution of differential equations using multiquadric radial basis function networks. Neutral Networks, 14: 185-199.
CrossRef  |  Direct Link  |  

23:  Aminataei, A. and M.M. Mazarei, 2005. Numerical solution of elliptic partial differential equarions using direct and indirect radial basis function networks. Eur. J. Sci. Res., 2: 5-5.

24:  Aminataei, A. and M.M. Mazarei, 2008. Numerical solution of Poisson's equation using radial basis function networks on the polar coordinate. Comput. Math. Appl., 11: 2887-2895.
CrossRef  |  Direct Link  |  

25:  Micchelli, C.A., 1986. Interpolation of scatter data: Distance matrices and conditionally positive definite functions. Constr. Approx., 2: 11-22.
CrossRef  |  Direct Link  |  

26:  Madych, W.R. and S.A. Nelson, 1990. Multivariable interpolation and conditionally positive definite functions. II. Math. Comput, 54: 211-230.
Direct Link  |  

27:  Kansa, E.J., 1990. Multiquadrics: A scattered data approximation scheme with applications to computational fluid-dynamics-I surface approximations and partial derivative estimates. Comput. Math. Applied, 19: 127-145.
CrossRef  |  Direct Link  |  

28:  Vanani, S.K. and A. Aminataei, 2009. Multiquadric approximation scheme on the numerical solution of delay differential systems of neutral type. Math. Comput. Modell., 49: 234-241.
CrossRef  |  

29:  Vanani, S.K. and A. Aminataei, 2008. On the numerical solution of neutral delay differential equations using multiquadric approximation scheme. Bull. Korean Math. Soc., 45: 663-670.
Direct Link  |  

30:  Vanani, S.K. and A. Aminataei, 2011. Numerical solution of differential algebraic equations using multiquadric approximation scheme. Math. Comput. Modell., 53: 659-666.
CrossRef  |  Direct Link  |  

31:  Madych, W.R., 1992. Miscellaneous error bounds for multiquadric and related interpolators. Comput. Math. Appl., 24: 121-138.
CrossRef  |  Direct Link  |  

32:  Arikoglu, A. and I. Ozkol, 2007. Solution of fractional differential equations by using differential transform method. Chaos Solitons Fractals, 34: 1473-1481.
CrossRef  |  Direct Link  |  

33:  Abdulaziz, O., I. Hashim and S. Momani, 2008. Application of homotopy perturbation method to fractional IVPs. J. Comput. Applied Math., 216: 574-584.
CrossRef  |  Direct Link  |  

34:  Kansa, E.J., 1990. Multiquadrics: A scattered data approximation scheme with applications to computational fluid-dynamics-II solutions to parabolic, hyperbolic and elliptic partial differential equations. Comput. Math. Appl., 19: 147-161.
CrossRef  |  

©  2021 Science Alert. All Rights Reserved