NUMERICAL SOLUTIONS OF VOLTERRA INTEGRO-DIFFERENTIAL EQUATIONS USING GENERAL LINEAR METHOD

. In this paper, a third order General Linear Method for ﬁnding the numerical solution of Volterra integro-diﬀerential equation is considered. The order conditions of the proposed method are derived based on techniques of B-series and rooted trees. The integral operator in Volterra integro-diﬀerential equation approximated using Simpsons rule and Lagrange interpolation is dis-cussed. To illustrate the eﬃciency of third order General Linear Method, we compare the method with a third order Runge-Kutta method.

1. Introduction. We consider solving the first order of Volterra integro-differential equation (VIDE) given in the form y (x) = F x, y(x), λ where f(x) and K(x,t) are known functions while y(x) is the solution to be obtained. Applications of VIDE arise in various fields such as biology, physics, engineering and more (see [11]). Consequently abundant numerical methods have appeared on finding the solutions to problem 1. Raftari in [10] solve the VIDE using Homotopy perturbation method along with trapezoidal rule for the integral operator and compared it with the finite difference method. Followed by Zarebnia in [12], the author considered the numerical solutions of VIDE by means of the Sinc collocation method where Sinc methods are direct solvers to integral equations. Later, Filiz in [5] applied several Runge-Kutta methods with different order associated with trapezoidal rule and Simpsons rule to VIDE. Then, Filiz in [6] used a higher order method which is the Runge-Kutta-Fehlberg method to solve the VIDE. He developed a new numerical routine for the integral part by using Lagrange interpolation and combination of various numerical quadrature rules to attain higher accuracy.
In this work, we focused on the General Linear Method (GLM) as a numerical method to obtain the solutions of VIDE. GLM was first introduced by Butcher in [1] as solver to initial value problems of ordinary differential equations. The GLM is known as natural generalizations of the classical Runge-Kutta and linear multistep methods. The combination of these methods is to achieve a high accuracy with high stage order method while incorporates good stability. Moreover, the order conditions of GLM are derived based on techniques of B-series and rooted trees (see [7], [3]). These techniques are alternatives to conventional way which is the Taylor method for finding order conditions to Runge-Kutta method and convenient for higher order method. Recently, Rabiei et al. in [9] proposed the General Linear Method for solving fuzzy differential equations and based on that, in this paper we construct the GLM adapted with suitable integration rules to solve VIDE problems.
The rest of the paper is organized as follows. In Section 2, we presented the general formula of GLM. In Section 3, implementation of GLM with suitable integration rules for Volterra integro-differential equation is given. Finally, numerical results and graphical illustrations are showed to test the performance of the proposed method.
2. General Linear Method. Consider the initial value problem of first order ordinary differential equation as follows: y (x) = f (x, y(x)), y(x 0 ) = y 0 (2) The general formulae for General Linear Method presented by Butcher in [1], [7] is given by: We denote the output approximations from step number n by y [n] i , i = 1, 2, , r, the stage values by Y i = 1, 2, ,s and the stage derivatives by F i = 1, 2, , s. Simply we have: Butcher's Tableau is used to represent the coefficients of GLM as given in Table 1.
The derivation of GLM can be further learned by understanding the concept of rooted trees. These rooted trees are first introduced by Butcher in [2] and used to obtain the order conditions which are needed to derive the GLM. Rabiei et al. in [9] use trees of order one to four. Aside from that, associated closely with the rooted trees is the composition law. This law is important in the construction of  [3]). The tedious calculations of Taylor series expansion can be avoided by conveniently apply this law in the derivation of algebraic order conditions. Definition 2.1. (see [1], [7]) Let a : T # → R by a mapping, then the form of a formal series is given as follows: is the B-series, where T # = denoting the set of rooted trees together with an additional empty tree ∅, a(∅) = 1, |t| = order of tree, γ(t) = density of tree, σ(t) = symmetry of tree, In order to derive the order conditions for GLM, we need to satisfy theorem below: [1]) The general linear method (A,U,B,V) has order p if there exists ξ ∈ X r and η ∈ X s 1 , such that, for every tree t satisfying |t| ≤ p, From Theorem 2.2, η represents the stages of method while Eξ is the output approximations computed after a time-step. The general coefficients of GLM of order three with r = 2, s = 3 are given in Table 2. Table 2. Matrix coefficients of GLM with s = 3, r = 2.
Based on Theorem 2.2, Rabiei et al. in [9] obtained the order conditions for a third order GLM as given in Table 3.  Table 3. Order conditions of GLM up to order three.
No Order conditions Based on the order conditions presented in Table 3, Rabiei et al. in [9] achieved the coefficients of third order with three stages of GLM method as shown in Table  4. In addition, we solve the order conditions in Table 3 to obtain another two new sets of coefficients of the third order GLM as presented in Tables 5and 6.
Based on the existing coefficients set and two new sets of GLM, we compute the error norms for all sets. The principal error norm is expressed as follows (see [4]).
The order conditions which are not satisfied by coefficients of GLM are as follows: Based on the error norms obtained in equations 9 -11, it is clear that the Set 1 of coefficients of GLM have smaller error norm compared with two other sets. Therefore, we use Set 1 coefficients for the following sections.
3. Implementation. Let us denote the integral operator in equation 1 as follows then the VIDE is given as The third order GLM is combined with appropriate numerical quadrature rule in order to approximate the integral part given in equation 12. To evaluate the integral z(x), we use the combination of composite Simpsons II rule (see [8]) for the interval [x 0 , x n ] and Lagrange interpolation for interval [x n , x ( n + c i )] at points x = −1, x = 0, x = c i . The third order General Linear Method applied to approximate equation 13 on a set of interval [0, X] of equally spaced grid points x 0 < x 1 .. < x N = X where 0 ≤ n ≤ N with step size h = (X−x0) N may be written as: For interval [x 0 , x n ], if n = 1, trapezoidal rule is used then else if n is even, composite 1/3 Simpsons rule is used then else if n is odd, composite Simpsons II rule is used then Therefore General Linear Method after employing composite Simpsons rule and Lagrange interpolation is given as follows hF (x n , Y 1 , z n ) + a 32 hF (x n+ 1 2 , Y 2 , z n+ 1 2 ) + u 31 y 1 (x n ) + u 32 y 2 (x n ), Step size s Number of stage MAXE Maximum error GLM Third order General Linear Method proposed in this paper RK Third order Runge-Kutta method by Butcher in [2] TFE Total function evaluations Problem (1). Consider the Volterra integro-differential equation (see Filiz [5]) y (x) = 1 − x 0 y(t)dt, y(0) = 0, with exact solution y(x) = sin?(x). Problem (2). Consider the Volterra integro-differential equation (see Filiz [5]) y (x) = 1 + x 0 y(t)dt, y(0) = 0, with exact solution y(x) = sinh?(x). Problem (3). Consider the Volterra integro-differential equation (see Zarebnia [12]) Problem (4). Consider the Volterra integro-differential equation (see Zarebnia [12]) y (x) = −y(x)+ x 0 e (t−x) y(t)dt, y(0) = 1, with exact solution y(x) = exp(−x)cosh(x). Problem (5). Consider the Volterra integro-differential equation (see Zarebnia [12]) x 0 x t+1 y(t)dt, y(0) = 0 with exact solution y(x) = ln(x + 1).  The result in Table 7 for problem 1 showed that the maximum error computed by the third order GLM with error exponent of 10 −6 is similar with the RK method. However, as the step size decreases we observed that the maximum errors of GLM are lower than RK method indicating that the accuracy of GLM for approximating the Volterra integro-differential equation is better than that of RK. In this case, the error of GLM came out as ×10 −13 which is better than the error of ×10 −12 by RK     Table 8 and 9 for problems 2 and 3 further prove that our proposed method clearly outperforms the RK method by having smaller errors. On top of that, the solutions obtained by GLM in Table 10 and 11 for problems 4 and 5 apparently outshine the RK method as well. Table 12 shows the total number of function evaluation used for GLM and RK GLM. GLM requires one step method to obtain the initial value for the second solution of GLM y 2 . We use fourth order Runge-Kutta method (RK4) to obtain that initial value. Therefore GLM should has extra function evaluation than RK. However in solving the VIDE, we use 3 point Lagrange interpolation for both GLM and RK. Therefore we use RK4 to obtain the initial value at y(x = −1) for both GLM and RK. Take note that for the GLM, we take y  . Conclusion. In this paper, we constructed the third order General Linear Method adapted with Simpsons II rule. Handling the integral operator in Volterra integro-differential equations using Lagrange interpolation is demonstrated as well. Then we applied the GLM on some test problems and compared the results with RK method. Numerical results showed that the GLM is more accurate than RK method in solving VIDE. The GLM showed the better performance compared to the same order RK method for all tested problems. It is mainly due to the structure of GLM that naturally generalizes the linear multistep method. It can be caused by the existence of the second output component, y