^{1}

^{2}

^{3}

Our focus is the development and implementation of a new two-step hybrid method for the direct solution of general second order ordinary differential equation. Power series is adopted as the basis function in the development of the method and the arising differential system of equations is collocated at all grid and off-grid points. The resulting equation is interpolated at selected points. We then analyzed the resulting scheme for its basic properties. Numerical examples were taken to illustrate the efficiency of the method. The results obtained converge closely with the exact solutions.

We consider the numerical solution of initial value problem of the form:

y ″ = f ( x , y , y ′ ) ; y ( a ) = y ( 0 ) ; y ′ ( a ) = γ (1)

In practice, higher order ordinary differential equations of this form y n = f ( x , y ′ , y ″ , ⋯ , y n − 1 ) , is solved by reducing it to systems of first order differential equation of the form:

y ′ = f ( x , y ) , y ( a ) = 0 , f ∈ c [ a , b ] , x , y ∈ ℝ n (2)

then an approximate method is applied to solve the resulting Equation (2) as widely discussed by Fatunla [

Although some of the aforementioned authors have made use of Taylor series, but little has been said with the use of Taylor series as a major method of implementation. So, Our idea is to use Taylor series algorithm to evaluate

y n + j , y ′ n + j , j = 1 , 2 and y n + v , y ′ n + v , v = 1 2 , 1 3 , 2 3 , 3 2 , ⋯ and calculate f ′ , f ″ by the

use partial derivative technique. Thus, two-step hybrid methods in the Taylor series mode are developed to solve second order ordinary differential equations directly.

In this section, power series is considered as an approximate solution to the general second order problems:

f ( x , y , y ′ , y ″ ) = 0 ; y ( a ) = y ( 0 ) ; y ′ ( a ) = γ (3)

of the form:

y ( x ) = ∑ j = 0 2 k + 1 a j x j (4)

The first and second derivative of (3) are respectively given as:

y ′ ( x ) = ∑ j = 1 2 k + 1 j a j x j − 1 (5)

y ″ = ∑ j = 2 2 k + 1 j ( j − 1 ) a j x j − 2 (6)

Combining (2) and (5), we generate the differential system

∑ j = 2 2 k + 1 j ( j − 1 ) a j x j − 2 = f ( x , y , y ′ ) , (7)

we develop the hybrid scheme using (3) and (5) as interpolation and collocation equations in this work.

Collocating (6) at selected grid and off-grid points, x = x n + 1 , 0 ≤ i ≤ 2 and interpolating (3) at selected grid and off-grid points, it results into a system of equations:

∑ j = 2 2 k + 1 j ( j − 1 ) a j x j − 2 = f n + i , 0 ≤ i ≤ 2 (8)

∑ j = 2 2 k + 1 a j x j = y n + i , 0 ≤ i ≤ 2 (9)

where, x n + i = x n + i h , solving Equations ((7) and (8)), a ′ j s , yield a method ex-

pressed in the form:

y k ( x ) = ∑ j = 0 k α j ( x ) y n + j + ∑ j = 0 k β j ( x ) f n + j , (10)

where k = 2 and f n + j = f ( x n + j , y n + j , y ′ n + j ) , 0 ≤ 2

It implies

a 0 + a 1 x n + a 2 x n 2 + a 3 x n 3 + a 4 x n 4 + a 5 x n 5 + a 6 x n 6 + a 7 x n 7 + a 8 x n 8 = y n (11)

a 0 + a 1 x n + 1 + a 2 x n + 1 2 + a 3 x n + 1 3 + a 4 x n + 1 4 + a 5 x n + 1 5 + a 6 x n + 1 6 + a 7 x n + 1 7 + a 8 x n + 1 8 = y n + 1 (12)

2 a 2 + 6 a 3 x n + 12 a 4 x n 2 + 20 a 5 x n 3 + 30 a 6 x n 4 + 42 a 7 x n 5 + 56 a 8 x n 6 = f n (13)

2 a 2 + 6 a 3 x n + 1 3 + 12 a 4 x n + 1 3 2 + 20 a 5 x n + 1 3 3 + 30 a 6 x n + 1 3 4 + 42 a 7 x n + 1 3 5 + 56 a 8 x n + 1 3 6 = f n + 1 3 (14)

2 a 2 + 6 a 3 x n + 2 3 + 12 a 4 x n + 2 3 2 + 20 a 5 x n + 2 3 3 + 30 a 6 x n + 2 3 4 + 42 a 7 x n + 2 3 5 + 56 a 8 x n + 2 3 6 = f n + 2 3 (15)

2 a 2 + 6 a 3 x n + 1 + 12 a 4 x n + 1 2 + 20 a 5 x n + 1 3 + 30 a 6 x n + 1 4 + 42 a 7 x n + 1 5 + 56 a 8 x n + 1 6 = f n + 1 (16)

2 a 2 + 6 a 3 x n + 4 3 + 12 a 4 x n + 4 3 2 + 20 a 5 x n + 4 3 3 + 30 a 6 x n + 4 3 4 + 42 a 7 x n + 4 3 5 + 56 a 8 x n + 4 3 6 = f n + 4 3 (17)

2 a 2 + 6 a 3 x n + 5 3 + 12 a 4 x n + 5 3 2 + 20 a 5 x n + 5 3 3 + 30 a 6 x n + 5 3 4 + 42 a 7 x n + 5 3 5 + 56 a 8 x n + 5 3 6 = f n + 5 3 (18)

2 a 2 + 6 a 3 x n + 2 + 12 a 4 x n + 2 2 + 20 a 5 x n + 2 3 + 30 a 6 x n + 2 4 + 42 a 7 x n + 2 5 + 56 a 8 x n + 2 5 = f n + 2 (19)

Writing these system of equations in matrix form:

Using Gaussian elimination method, the unknown coefficients a ′ j s can be obtained. Putting a ′ j s back into (3) gives (10):

The coefficients α i ′ s ( t ) , β j ′ s ( t ) are continuous coefficients obtained using the transformation t = 1 h ( x − x n + k − 1 ) , t ∈ ( 0,1 ]

d t d x = 1 h .

Then simplifying the continuous α j ′ s , β j ′ s , and taking their first derivatives, we have:

α 0 ( t ) ′ = − 1 h ,

α 1 ( t ) ′ = − 1 h ,

β 0 ( t ) ′ = 47 h 13440 ,

β 1 3 ( t ) ′ = 327 h 2240 ,

β 2 3 ( t ) ′ = 111 h 890 ,

β 1 ( t ) ′ = 1088 h 3360 ,

β 4 3 ( t ) ′ = 93 h 640 ,

β 5 3 ( t ) ′ = 1095 h 2240 ,

β 2 ( t ) ′ = 1359 h 13440 .

Then, putting t = 1 gives:

y n + 2 = 2 y n + 1 + y n + h 2 6720 { 47 f n + 2 + 810 f n + 5 3 + 1377 f n + 4 3 + 2252 f n + 1 + 1377 f n + 2 3 + 810 f n + 1 3 + 857 f n } (21)

its first derivative

y ′ n + 2 = 1 h [ y n + 1 − y n ] + h 2 6720 { 1359 f n + 2 + 6570 f n + 5 3 + 1953 f n + 4 3 + 4352 f n + 1 + 1665 f n + 2 3 + 1962 f n + 1 3 + 47 f n } (22)

with the order p = 8 , error constant C 10 = − 0.0069941 , and interval of absolute stability X ( Θ ) = ( − 14.1608 , 0 ) Implementation of the method using Taylor series algorithm to evaluate

y n + j , y ′ n + j , y n + v , y ′ n + v , f n + v , f n + j ,

where, j ′ s = 1 , 2 and v ′ s = 1 3 , 2 3 , 4 3 , 5 3 and,

f n + v = f ( x n + v , y n + v , y ′ n + v ) ,

such that

y n + v = y n + v h y ′ n + ( v h ) 2 2 ! f n + ( v h ) 3 3 ! f ′ n + ( v h ) 4 4 ! f ″ n + ⋯ (23)

and,

y ′ n + v = y ′ n + v h f n + ( v h ) 2 2 ! f ′ n + ( v h ) 3 3 ! f ″ n + ( v h ) 4 4 ! f ‴ n + ⋯ (24)

Also,

f n + j = y ″ ( x n + j h ) = f n + j h f ′ n + ( j h ) 2 2 ! f ″ n + ⋯ (25)

f n = f ( x n , y n , y ′ n )

f ( i ) = f ( i ) ( x n , y n , y ′ n ) , i = 1 , 2

Finding the partial derivative f ′ , f ″ , ⋯ as follows

d f d x = f ′ = δ f δ x + δ f δ x y ′ + δ f δ y ′ f (26)

f ″ = d 2 f d x 2 = 2 ( A y ′ + B f ) + C f y ′ + D + E , (27)

where,

A = ∂ 2 f ∂ x ∂ y + f ∂ 2 f ∂ y ∂ y ′ (28)

B = ∂ 2 f ∂ x ∂ y ′ (29)

C = ∂ f ∂ x + y ′ ∂ f ∂ y + f ∂ f ∂ y ′ (30)

D = ∂ 2 f ∂ x 2 + ( y ′ ) 2 ∂ 2 f ∂ y 2 + f 2 ∂ 2 f ∂ ( y ′ ) 2 (31)

E = f ∂ f ∂ y (32)

We shall consider the analysis of the basic properties of our methods which includes the order, the region of absolute stability and the zero stability of the methods.

The local truncation error with k-step linear multistep m method which is in line with Lambert (1973), is taken to be linear difference operator l defined by

l [ y ( x ) ; h ] = ∑ j = 0 k [ α j y ( x n + j ) − h β j y ( x n + j ) ] (33)

Thus, expanding (21) as Taylor series about point x and comparing coefficients of h k , the scheme will be of order p = 8 with error constant C p + 2 = − 0.0069941

L [ y ( x ) , h ] = C 0 y ( x n ) + C 1 y ′ ( x n ) + C 2 y ″ ( x n ) + ⋯ + C p y p ( x n ) , (34)

where C p , p = 0 , 1 , ⋯ are the constant coefficients given as:

C 0 = ∑ j = 0 k α j C 1 = ∑ j = 0 k j α j and C p = 1 p ! [ ∑ j = 0 k j α j − p ( p − 1 ) ( ∑ j = 0 k j p − 1 β j + ∑ j = 0 k q p − 1 β q j ) ] } (35)

In line with [

A linear multistep method is consistent if the following conditions are satisfied:

1) The order p ≥ 1 .

2) p ( 1 ) = 0 , p ′ ( 1 ) = σ ( 1 ) .

3) ∑ j = 0 k α j = 0 .

4) ∑ j = 0 k j α j = ∑ j = 0 k β j .

Equation (21) has its first characteristic polynomial to be:

ρ ( r ) = r 2 − 2 r + 1 (36)

The method is zero stable since they have roots r = 1 twice.

In order to establish the region of absolute stability, we apply the boundary locus method as in [

θ = ρ ( r ) δ ( r )

where,

r = e i θ = cos ( θ ) + i sin ( θ )

From scheme (21), we have: ρ ( r ) = r 2 − 2 r + 1

and

σ ( r ) = 1 6720 [ 47 r 2 + 810 r 5 3 + 1377 r 4 3 + 2252 r + 1377 r 2 3 + 810 r 1 3 + 857 ]

so that

h ( θ ) = ρ ( e i θ ) δ ( e i θ )

which implies

h ( θ ) = 1 6720 [ 47 r 2 + 810 r 5 3 + 1377 r 4 3 + 2252 r + 1377 r 2 3 + 810 r 1 3 + 857 ] (37)

h ( θ ) = [ cos ( 2 θ ) + i sin ( 2 θ ) − 2 cos ( θ ) − 2 i sin ( θ ) + 1 ] × 6720 ( 47 cos ( 2 θ ) + 47 i sin ( 2 θ ) + 810 cos ( 5 θ 3 ) + 810 i sin ( 5 θ 3 ) + 1377 cos ( 4 θ 3 ) + 1377 i sin ( 4 θ 3 ) + 2252 cos ( θ ) + 2252 i sin ( θ ) + 1377 cos ( 2 θ 3 ) + 1377 i sin ( 2 θ 3 ) + 810 cos ( 2 θ 3 ) + 810 i sin ( 2 θ 3 ) + 47 ) − 1

Considering the values of θ for 0 ≤ θ ≤ 180 at intervals of 30 θ gives the region of absolute stability to be ( − 14.1608,0 ) .

We test the accuracy of the proposed scheme on some numerical problems, and the results are compared with existing methods.

Problem 1:

y ″ = x ( y ′ ) 2 , y ( 0 ) = 1 , y ′ ( 0 ) = 0.5 , h = 0.1 32 (38)

Exact solution

y ( x ) = 1 + 1 2 log 10 ( 2 + x 2 − x )

The numerical results of the problem is shown in

Problem 2:

y ″ = ( − 6 / x ) y ′ − ( 4 / ( x ) 2 ) y , y ( 1 ) = 1 , y ′ ( 1 ) = 1 , h = 120 (39)

(x) | YEX | YC | ERRNew |
---|---|---|---|

0.2 | 1.100335347731075300 | 1.100335347731045300 | |

0.4 | 1.20273255405481600 | 1.11273255405480200 | |

0.6 | 1.309519604203111900 | 1.009519604203101000 | |

0.8 | 1.423648930193603500 | 1.123648930123598200 | |

1.0 | 1.549306144334058600 | 1.129306144334043400 |

(x) | YEX | YC | ERRNew |
---|---|---|---|

1.0125 | 1.0117410181167988400 | 1.011741018167989300 | |

1.0188 | 1.017066494235672900 | 1.017066494235672900 | |

1.0250 | 1.017066494235672900 | 1.022049163629432000 | |

0.8 | −1.2255409228492467900 | −1.225540922161721500 | |

1.0313 | 1.026703577500806200 | 1.026703577500806700 |

Note: YEX = Yexact, YC = Ycomputed, ERRNew = Error in new method.

Exact solution

y ( x ) = 1 − e x

The numerical results of the problem is shown in

A Linear Multistep method which implements a Taylor’s series algorithm is developed for the direct solution of general second order initial value problems of

ordinary differential equations without reduction to systems of first order differential equation. The derivatives of continuous scheme to any order were computed implementing Taylor’s series algorithm. The accuracy of the method was tested with two test problems, and results were compared with Awoyemi and Kayode [

Owolanke, A.O., Uwaheren, O. and Obarhua, F.O. (2017) An Eight Order Two-Step Taylor Series Algori- thm for the Numerical Solutions of Initial Value Problems of Second Order Ordinary Differential Equations. Open Access Lib- rary Journal, 4: e3486. https://doi.org/10.4236/oalib.1103486