An Order Four Continuous Numerical Method for Solving General Second Order Ordinary Differential Equations

Continuous hybrid methods are now recognized as efficient numerical methods for problems whose solutions have finite domains or cannot be solved analytically. In this work, the continuous hybrid numerical method for the solution of general second order initial value problems of ordinary differential equations is considered. The method of collocation of the differential system arising from the approximate solution to the problem is adopted using the power series as a basis function. The method is zero stable, consistent, convergent. It is suitable for both non-stiff and mildly-stiff problems and results were found to compete favorably with the existing methods in terms of accuracy. DOI:10.46481/jnsps.2021.150


Introduction
We consider the second order Ordinary Differential Equation (ODEs) y = f (x, y, y ), y(µ) = ω 0 , y (µ) = ω 1 (1) Equation (1) occur virtually every areas of physical or biological process in connection with numerous problems that are encountered in various aspects of everyday life. It is well conceived that this type of equation can either be solved directly or solved by reducing to system of first order differential equations before applying different methods available to solve the resulting system of first order ODEs Chan et al. [1], Gholamtabar and Parandin [2].
Various linear multistep methods with different order of accuracy have been developed for the solution of 1 which varies from discrete linear multistep methods to continuous linear multistep methods. Lambert and Watsan, [3] reported that linear multistep methods generally are more efficient in terms of accuracy with weak stability properties for a given number of evaluation per step and suffered the disadvantage of requiring additional starting values and special procedures for changing step length h. It is also good to note that continuous linear multistep methods have advantages over the discrete methods in such a way that they give better error estimation, they provide a simplified form of coefficients for further evaluation at different points, and they provide solutions at all interior points within the interval of integration Saravi and Mirrajei, [4], Kayode and Awoyemi, [5], Golbabai and Arabshahi [6]. Among the first methods developed are first derivative methods that are implemented in predictor-corrector mode, and Taylor series expan-Obarhua Adegboro / J. Nig. Soc. Phys. Sci. 3 (2021) 42-47 43 sion are adopted to provide the starting values. The identified setbacks of the predictor-corrector methods are; they are very costly to implement and reduced order of accuracy of the predictors. Recently, authors have proposed different methods of higher order differential equations to improve on the existing setbacks. Such improved methods are Kayode and Adeyeye, [7,8] and Kayode and Obarhua,[9,10]. They independently proposed linear multistep methods of higher order of accuracies and same order of main predictors and the correctors and hence improved significantly the accuracies of the methods.
This work proposed an accurate continuous numerical hybrid method for direct solution of initial value problems of ODEs. The derived method is capable to handle stiff, mildly stiff, nonlinear and engineering problems modeled as a second order initial value problem of ODEs.

Derivation of the Method
We define the general power series approximation solution in the form Equating (4) with (1) gives Equation (2) is interpolated at x n+i , i = 2, 5 2 and (5) is collocated at x n+c , c = 0(1)3. Therefore, interpolation and collocation equation at the selected grid and offstep points give rise to system of equations which can be express in matrix form Gaussian elimination method is then applied to solve equation (6) to obtain the unknown coefficients a j s which is then substituted into (2). Continuous system is obtained after some algebraic simplifications. [11], the continuous coefficients are obtained as follows The first derivatives of equation (7) are Evaluating equation (7) and (8) at t = 1 yield the discrete order continuous numerical scheme its first derivative is given as The values of r is taken in the interval r ∈ (2, 3) to obtain a particular discrete hybrid method. For the purpose of testing the properties of equation (9), the value of r is taken to 5 2 to have with its first derivative given by 3. Implementation of the Method (11) In order to implement the implicit linear one-point discrete scheme (11) and its derivative (12), the symmetric explicit schemes and their derivatives are also developed by the same procedure for the evaluation of y n+3 and y n+3 contained in f n+3 in the main scheme (11) and its derivative (12). (13) and its first derivative as Other explicit schemes were also generated to evaluate other starting values using Taylor series expansions to evaluate the values for y n+i , y n+i as y n+i = y n +( jh)y n + ( jh) 2 2! f n + ( jh) 3 3! ∂ f n ∂x n + y n ∂ f n ∂y n + f n ∂ f n ∂y n +o(h 4 ) (15) and y n+i = y n +( jh) f n + ( jh) 2 2! ∂ f n ∂x n + y n ∂ f n ∂y n + f n ∂ f n ∂y n +o(h 4 )(16)

Region of Absolute Stability
In other to investigate the periodic stability properties of the numerical methods for solving the initial-value problem equa-tion 1 and the interval of periodicity, Lambert and Watson [3] introduced the following scalar test problem as Based on the theory developed in Lambert and Watson [3], when multistep method is applied to the scalar test equation (17), a difference equation of the form is obtained, where H = ph, h is the step length and y n is the computed approximation to y(x 0 + nh), n = 0, 1, 2, . . . Then, we have following definitions.

Definition. (See Konguetsof and Simos, [12]) Numerical method (19) has an interval of periodicity
Definition. Following [3], a numerical method is P-stable if its interval of periodicity is (0, ∞). Therefore, we obtain the interval of periodicity of the new method, which is Equal to (0, -2.4) and the stability domain of the method is as shown in Figure 1.

Order and Error Constant of the Method
The method proposed by Lambert (1973) in Olanegan et al. [13] is adopted in this paper, with linear operator: We associate the linear operator L with the scheme and defined as Where α 0 and β 0 are both non-zero and y(x) is an arbitrary function which is continuous and differentiable on the interval [a, b]. Expanding the form y(x) and y (x) in Taylor series and comparing coefficients of h, we obtained ∆[y(x); h] = C 0 y(x) + C 1 hy (x) + · · · + C p h p y p (x) The method (11) and its associate linear difference operator (13) are said to have order p if c 0 = c 1 = ··· = c p+1 = 0 and c p+2 0. The value c p+2 is called error constant. Therefore, in this paper, the method (11)

Consistency of the Method
A numerical method is said to be consistent if the following conditions are satisfied 1. The order of the method must be greater or equal to 1, p ≥ 1.
Where ρ(r) and σ(r) are the first and second characteristics polynomial of our method. according to Adesanya et al. [14], the first condition is a sufficient condition for a method to be consistent. Since our method is of order 4 then it is consistent.

Convergence
A method is said to be convergent if and only if it is consistent and zero stable, hence our method is convergent.

Numerical Examples
The method is applied to solve the following linear and nonlinear second order initial value problems of ordinary differential equations directly without reduction to system of first order equations. Problem 1: y = y , y(0) = 0; y (0) = −1; h = 0.1 Theoretical solution: y(x) = 1 − e x This problem has been used in Kayode and Adeyeye [8] of order six to check the behavior of the methods. Table 1 shows the absolute errors of the methods for the purpose of comparison.
The obtained results in the Table give the good performance of the proposed method. Problem 2: y + 1001y + 1000y = 0, y(0) = 1; y (0) = −1; h = 0.05 Theoretical solution: y(x) = e −x The Problem 2 was solved by Anake [15] of order 4. The new method was applied to solve it for the purpose of comparison. The results are shown in Table 2. Problem 3: y = 100y , y(0) = 1; y (0) = −10; h = 0.01 Theoretical solution: Table shows the absolute errors at different points of the integration interval when h = 0.01was solved by Awari [16] of order five. The new method was applied to solve it for the purpose of comparison. The results show that the proposed method performed better than Awari [16].  Table 5.