A FAMILY OF EXPONENTIALLY FITTED MULTIDERIVATIVE METHOD FOR STIFF DIFFERENTIAL EQUATIONS

In this paper, an A-stable exponentially fitted predictor-corrector using multiderivative linear multistep method for solving stiff differential equations is developed. The method which is a two-step third derivative method of order five contains free parameters. The numerical stability analysis of the method was discussed, and found to be A-stable. Numerical examples are provided to show the efficiency of the method when compared with existing methods in the literature that have solved the set of problems.


INTRODUCTION
In order to solve stiff initial value problems in ordinary differential equations efficiently, many new methods have been developed in recent years, which satisfy certain stability requirements.
The property of A-stability is desirable in formulas to be used in the solution of stiff systems of differential equations especially from chemical kinetics and the discretisation of partial differential equations. Dahlquist [7], proved that A-stable linear multistep formulas must be implicit, its maximum order is two and of those of second order, the one with the smallest truncation error coefficient is the trapezoidal rule. In nearly every linear systems of differential equation, which have widely disperse eignevalues, high order A-stable formulas are particularly appropriate since they allow integration to proceed with a larger step size.
Thus, the need to develop a high order A-stable implicit multistep formula, which uses linear combinations of derivatives higher than the first, give rise to the development of multiderivative multistep formulas. Lambert [11], and Enright [8], pointed out that Multiderivative methods give high accuracy and possess good stability properties when used to solve first order initial value problems in ordinary differential equation.

THE GENERAL MULTIDERIVATIVE MULTISTEP METHOD
The general multiderivative multistep method is given by,

∑ ∑ ∑
Where is the derivative of evaluated at , and are real constants with and is the appropriate numerical solution evaluated at the point . In order to remove the arbitrary constant in (2.1) we shall always assume that , and ∑ | | ∑ | |

DEVELOPMENT OF NEW EXPONENTIALLY FITTED MULTIDERIVATIVE METHODS
A numerical integration formula is said to be exponentially fitted at (complex) value if when the method is applied to the scalar test problem , with exact initial condition, the characteristic equation satisfies the relation .
However, the idea of using exponentially fitted formulas for the appropriate numerical integration of certain classes of stiff systems of first order ordinary differential equation of the form, which was originally proposed by Liniger and Willoughby [12], is to derive integration formulas containing free parameters (other than the step length of integration) and then to choose these parameters so that a given function , where is real, satisfies the integration formula exactly.
Liniger and Willoughby [12] derived three 1-step integration formulas with orders ranging from 1 -3. Their results revealed that for all choices of the fitting parameter , their formulas are A-stable.

DERIVATION OF FAMILY OF TWO STEP EXPONENTIALLY FITTED METHODS
The objective of this paper is to develop a two-step, third derivative multiderivative exponentially fitted formulas. (i.e. and ).
For this purpose, equation ( are respectively the first, second and third derivatives of . When we are deriving exponentially fitted multistep methods, the approach is to allow both (2.3) and (2.4) to posses free parameters order than the mesh size which allow it to be fitted automatically to exponential function.

DERIVATION OF METHOD OF FIFTH ORDER FORMULAS
The derivation of predictor-corrector integrations formula of order 5 involves two stages as was done in higher order. First we derived the order four predictor by setting in equation (2.3) to obtain five set of equations with 12 unknown parameters and we then obtain the following set of simultaneous equations When we solve the above equations we obtain, When these values of the parameters are substituted into (2.3), we obtain the predictor formula as: Now, for exponential fitting purpose, we apply (2.7) to scalar test function to obtain equation (2.8) below. Again to obtain the corresponding order 5 corrector formula, we obtain six set of simultaneous equation from (2.4)asfollows, 3 We impose the same condition as in predictor, and in addition, we let as free parameters , the values of the unknown parameters are obtained from (2.11) as, When these values are substitute into (2.4), we obtain the fifth order corrector formula as,

STABILITY CONSDIERATION OF THE METHOD
To examine the stability conditions required by this method, it is expected by maximum modulus theorem that the stability function of the method given by (2.16) satisfies| | .
In order to determine the interval of absolute stability of the method, we find limits of both and as and Similarly, from (2.16) we have, Thus we formed that ( ) and ( ) Now, we further verify analytically that the ranges of value of a andb represent the region of absolute stability of the new method. We illustrate this by taking a large sample S as shown in table 3.1 below;  However, for purpose of comparative analysis on the performance of the new scheme we denote AL 5 as the new method, CH4, CH5-Cash [6] method of order 4 and 5 respectively, J-K Jackson and Kenue [10], OK6-Okunuga [13], AB7, AB8, NM9 represent Abhulimen and Otunta [3] method of order 7, 8 and 9 respectively. F 5 Abhulimen[4] method of order 5, AF5-Abhulimen and Okunuga [1]and AG6-Abhulimen and Omeike [2].

Example 1
Non-linear stiff problems (Enright and Pryce) [9]  Step length h Method The system (4.2) can be rewritten as a firstorder system; Thus we obtain a 2 x 2 system of stiff IVP. The eigenvalues of the Jacobian matrix The general solution of (4.2) is . If we impose the initial conditions in the exact solution is .
The result of this problem using the newly derived methods are obtained at ; as given in table 4.2 below.
The exact solution is given as, For comparison purpose , we have the following as; AB7, AB8 and NM9 to represent Abhulimen and Otunta " Two step third derivative methods order seven, eight and nine" respectively. F 5 denote three-step second derivative scheme. As shown in the Table 3.1 above, the proposed method in this paper perform better than existingmethods in terms of accuracy.

CONCLUSION
The aim of this paper was to develop numerical method which provides solution to initial value problems with stiff differential equations via exponentially fitted integrators.Numerical experiments have been carried out using appropriate step size as required by each problem. Such problems which are stiff require small step size before the solution can be smooth. In general, the results from numerical experiment so presented in this paper, show that the new method perform effectively well when compared similar methods in the literature. Hence the aim and objective of this paper have been achieved.