A robust diagonally implicit block method for solving first order stiff IVP of ODEs

22  Download (0)

Full text

(1)

A Robust Diagonally Implicit Block Method for Solving First Order Stiff IVP of ODEs

Abdu Masanawa Sagir1 and Muhammad Abdullahi1*

1Department of Mathematical Sciences, Federal University Dutsin-Ma, Katsina State, Nigeria.

*Corresponding author: maunwala@gmail.com

Received: 13 September 2022; Accepted: 21 October 2022; Available online (In press): 31 October 2022

ABSTRACT

In this work, a block of diagonally implicit backward differentiation method with two off-step points for solving first order stiff initial value problem of ordinary differential equation was derived. In the proposed block method two approximate solution values of ๐‘ฆ๐‘›+1 and ๐‘ฆ๐‘›+2 with two off-step points ๐‘ฆ๐‘›+1

2 and ๐‘ฆ๐‘›+3

2 are computed concurrently for each iteration. The properties of the newly proposed method were found to be an A-stable, Zero stable and capable for solving first order Stiff IVPs. To validate the performance of the proposed method, some first order stiff IVPs are solved and the result obtained was compared with other existing numerical schemes.

From the tabulated results and the graphs plotted, the proposed method has shown advantages of accuracy in the scale error over the three methods and an advantage of executional time over two of the existing methods considered.

Keywords: Block, Diagonally, Implicit, Ordinary Differential Equation, Stiff IVPs

1 INTRODUCTION

In physical, social and life sciences so many real life problems can be modelled into equation, most often differential equation. Problems in electrical circuits, chemical reactions, mechanics, vibrations, and kinetics and population growth all can be modelled as differential equations and categorized into stiff and non-stiff. A stiff equation is a differential equation for which certain numerical methods for solving the equation are numerically unstable, unless the step size is taken to be extremely small. So, with regard to modern scientific and engineering cases most of the modelled equations derived turns to be stiff problems. Researchers are formulating various methods to obtain analytic and/or numerical solutions of the modelled stiff IVPs. But, stiff problem usually deviates from been solved analytically due its complexities and other phenomena which is found within its solution, the transient and steady state components found in its solution make explicit method difficult to handle with appreciated results. While, numerical solution is much more easier and obtainable in any form stiff IVP of ODEs. Most of the stiff cases have no analytical solutions at all. Hence, preferences are always channels to numerical methods that would solve any sort of stiff IVP of ODEs. The ultimate goal is to get a method with a solution that has absolutely minimum scale error and computational time.

(2)

The early work integration of stiff equations [1], the extended backward differential formula and modified extended BDF [2, 3]. Block method was given more recognition through the work of [4, 8, 9, 10, 11, 12]. New fifth order implicit block method for solving first order stiff ordinary differential equations, a new super class of block backward differentiation formula for stiff ordinary differential equations, an accurate block solver for stiff initial value problems [5, 6, 7]. More numerical solutions and Hybrid multistep block methods were developed [19, 20, 21, 22, 23, 28]. Numerical treatment of block method for the solution of ODEs, on the approximate solution of continuous coefficients for solving third order ordinary differential equations, an accurate computation of block hybrid method for solving stiff ordinary differential equations [13, 14, 15]. An A-stable block integrator scheme for the solution of first order system of IVPs of ordinary differential equations, order and convergence of the enhanced 3-point fully implicit super class of block backward differentiation formula for solving first order stiff initial value problems [16, 17]; extended 3-Point super class of block backward differentiation formula for solving first order stiff initial value problems [18]. Other works include the following [24], [25], [26], [27], [29], [30], [31], [32], [33], [34] and [35] all the methods highlighted above demonstrates very good stability properties, at one point or the other, with appreciated results in terms of accuracy and computational time.

This paper considers derivation and implementation of a robust diagonally implicit block method with two off-step points for solving a system of first order initial value problem of ordinary differential equations. The proposed method will generate newly equally spaced solution values concurrently. In addition, the method uses of a lower triangular matrix with identical diagonal entries, as such the coefficients of the upper triangular matrix entries are zero. The method is of the form

๐‘ฆโ€ฒ = ๐‘“(๐‘ฅ, ๐‘Œฬ‚), ๐‘Œฬ‚(๐‘Ž) = ๐œ‘ีฒ, ๐‘Ž โ‰ค ๐‘ฅ โ‰ค ๐‘ (1)

where ๐‘Œฬ‚ = (๐‘ฆ1, ๐‘ฆ2, ๐‘ฆ3, โ€ฆ โ€ฆ โ€ฆ ๐‘ฆ๐‘›), ีฒ๐œ‘ฬ… = (๐œ‘ีฒ1, ๐œ‘ีฒ2, ๐œ‘ีฒ3, โ€ฆ , ๐œ‘ีฒ๐‘›).

2 DERIVATION OF A ROBUST 2-POINT BBDF WITH OFF-STEP POINTS

In this section, two approximate solution values ๐‘ฆ๐‘›+1 and ๐‘ฆ๐‘›+2 with step size โ„Ž, and two off-step points ๐‘ฆ๐‘›+1

2

and ๐‘ฆ๐‘›+3 2

which are chosen at the points where the step size is halved are formulated in a block simultaneously. The formulae are computed using two back values ๐‘ฆ๐‘› and ๐‘ฆ๐‘›โˆ’1 with step size h. The formulae are derived with the aid of a diagram as shown in Figure 1.

(3)

Figure 1: Diagram for RDIBM derivation.

The proposed method (RDIBM) is of the form

โˆ‘1+๐‘˜๐‘—=๐‘œ๐›ผ๐‘—,๐‘–๐‘ฆ๐‘›+๐‘—โˆ’1= โ„Ž๐›ฝ๐‘˜,๐‘–[๐‘“๐‘›+๐‘˜+ ๐‘“๐‘›+๐‘˜โˆ’1]๐‘˜ = ๐‘– =1

2, 1,3

2, 2 (2)

where ๐‘˜ and ๐‘– have the same value. The formula (2) is derived using Taylorโ€™s series expansion about ๐‘ฅ๐‘›

Definition 2.1. According to [26], the linear operator ๐ฟ๐‘– associated with first, second, third and fourth point of the RDIBM with off-step points method is defined as follows:

๐ฟ๏ฎ[๐‘ฆ(๐‘ฅ๐‘›), โ„Ž]: ๐›ผ0,๏ฎ ๐‘ฆ๐‘›โˆ’1+ ๐›ผ1,๏ฎ ๐‘ฆ๐‘›+ ๐›ผ3

2,๏ฎ๐‘ฆ๐‘›+1

2

โˆ’ โ„Ž๐›ฝ๏น,๏ฎ [๐‘“๐‘›+๏น + ๐‘“๐‘›+๏นโˆ’1] = 0

๐ฟ๏ฎ [๐‘ฆ(๐‘ฅ๐‘›), โ„Ž]: ๐›ผ0,๏ฎ ๐‘ฆ๐‘›โˆ’1+ ๐›ผ1,๏ฎ ๐‘ฆ๐‘›+ ๐›ผ3

2,๏ฎ๐‘ฆ๐‘›+1

2

+ ๐›ผ2,๏ฎ ๐‘ฆ๐‘›+1โˆ’ โ„Ž๐›ฝ๏น,๏ฎ [ ๐‘“๐‘›+๏น

+๐‘“๐‘›+๏น โˆ’1] = 0

๐ฟ๏ฎ [๐‘ฆ(๐‘ฅ๐‘›), โ„Ž]: ๐›ผ0,๏ฎ๐‘ฆ๐‘›โˆ’1+ ๐›ผ1,๏ฎ ๐‘ฆ๐‘›+ ๐›ผ3

2,๏ฎ ๐‘ฆ๐‘›+1

2

+ ๐›ผ2,๏ฎ๐‘ฆ๐‘›+1+ ๐›ผ5

2,๏ฎ ๐‘ฆ๐‘›+3

2

โˆ’ โ„Ž๐›ฝ๏น,๏ฎ [ ๐‘“๐‘›+๏น

+๐‘“๐‘›+๏น โˆ’1] = 0 ๐ฟ๏ฎ[๐‘ฆ(๐‘ฅ๐‘›), โ„Ž]: ๐›ผ0,๏ฎ ๐‘ฆ๐‘›โˆ’1+ ๐›ผ1,๏ฎ ๐‘ฆ๐‘›+ ๐›ผ3

2,๏ฎ๐‘ฆ๐‘›+1

2

+ ๐›ผ2,๏ฎ ๐‘ฆ๐‘›+1+ ๐›ผ5

2,๏ฎ๐‘ฆ๐‘›+3

2

+ ๐›ผ3,๏ฎ ๐‘ฆ๐‘›+2

โˆ’โ„Ž๐›ฝ๏น ,๏ฎ [๐‘“๐‘›+๏น + ๐‘“๐‘›+๏น โˆ’1] = 0 }

(3)

Consider the following value of ๏น ๏ฎ& 'svalue in (3) for the cases below:

For cases 1, 2, 3 and 4 as in ๏น =๏ฎ =1

2,๏น = ๏ฎ = 1,๏น =๏ฎ =3

2 & ๏น =๏ฎ = 2 for the first, second, third and fourth point respectively, with the associated operator (๐ฟ1

2

, ๐ฟ1 , ๐ฟ3 2

& ๐ฟ2) related to (3) is written as

(4)

๐›ผ0,1

2

๐‘ฆ(๐‘ฅ๐‘›โˆ’ โ„Ž) + ๐›ผ1,1

2

๐‘ฆ(๐‘ฅ๐‘›) + ๐›ผ3

2,1

2

๐‘ฆ(๐‘ฅ๐‘›+1

2โ„Ž) โˆ’ โ„Ž๐›ฝ1

2,1

2

[๐‘“ (๐‘ฅ๐‘›+1

2โ„Ž) + ๐‘“ (๐‘ฅ๐‘›โˆ’1

2โ„Ž)] = 0 ๐›ผ0,1๐‘ฆ(๐‘ฅ๐‘›โˆ’ โ„Ž) + ๐›ผ1,1๐‘ฆ(๐‘ฅ๐‘›) + ๐›ผ3

2,1๐‘ฆ (๐‘ฅ๐‘›+1

2โ„Ž) + ๐›ผ2,1๐‘ฆ(๐‘ฅ๐‘›+ โ„Ž)

โˆ’โ„Ž๐›ฝ1,1[๐‘“(๐‘ฅ๐‘›+ โ„Ž) + ๐‘“(๐‘ฅ๐‘›)] = 0 ๐›ผ0,3

2

๐‘ฆ(๐‘ฅ๐‘›โˆ’ โ„Ž) + ๐›ผ13

2,๐‘ฆ(๐‘ฅ๐‘›) + ๐›ผ3 2,3

2

๐‘ฆ (๐‘ฅ๐‘›+1

2โ„Ž) + ๐›ผ2,3 2

๐‘ฆ(๐‘ฅ๐‘›+ โ„Ž) + ๐›ผ5 2,3

2

๐‘ฆ (๐‘ฅ๐‘›+3

2โ„Ž)

โˆ’โ„Ž๐›ฝ3

2,3

2

[๐‘“ (๐‘ฅ๐‘›+3

2โ„Ž) + ๐‘“ (๐‘ฅ๐‘›+1

2โ„Ž)] = 0 ๐›ผ0,2๐‘ฆ(๐‘ฅ๐‘›โˆ’ โ„Ž) + ๐›ผ1,2๐‘ฆ(๐‘ฅ๐‘›) + ๐›ผ3

2,2๐‘ฆ (๐‘ฅ๐‘›+12โ„Ž) + ๐›ผ2,2๐‘ฆ(๐‘ฅ๐‘›+ โ„Ž) + ๐›ผ5

2,2๐‘ฆ (๐‘ฅ๐‘›+32โ„Ž) +๐›ผ3,2๐‘ฆ(๐‘ฅ๐‘›+ 2โ„Ž) โˆ’ โ„Ž๐›ฝ2,2[๐‘“(๐‘ฅ๐‘›+ 2โ„Ž) + ๐‘“(๐‘ฅ๐‘›+ โ„Ž)] = 0 }

(4)

Expanding (๐‘ฅ๐‘›โˆ’ โ„Ž) , ๐‘ฆ(๐‘ฅ๐‘›) , ๐‘ฆ(๐‘ฅ๐‘›+1

2โ„Ž),๐‘ฆ(๐‘ฅ๐‘›+ โ„Ž), ๐‘ฆ(๐‘ฅ๐‘›+3

2โ„Ž), ๐‘ฆ(๐‘ฅ๐‘›+ 2โ„Ž),๐‘“ (๐‘ฅ๐‘›+1

2โ„Ž) , ๐‘“ (๐‘ฅ๐‘›โˆ’

1

2โ„Ž), ๐‘“ (๐‘ฅ๐‘›+32โ„Ž) , ๐‘“(๐‘ฅ๐‘›+ โ„Ž), ๐‘“(๐‘ฅ๐‘›+ 2โ„Ž) in (4) with a Taylorโ€™s series expansion about ๐‘ฅ๐‘› and collect the like terms gives

๐ถ0,1

2

๐‘ฆ(๐‘ฅ๐‘›) + ๐ถ1,1 2

โ„Ž๐‘ฆโ€ฒ(๐‘ฅ๐‘›) + ๐ถ3 2,1

2

โ„Ž2๐‘ฆโ€ฒโ€ฒ(๐‘ฅ๐‘›) + โ‹ฏ = 0 ๐ถ0,1๐‘ฆ(๐‘ฅ๐‘›) + ๐ถ1,1โ„Ž๐‘ฆโ€ฒ(๐‘ฅ๐‘›) + ๐ถ3

2,1โ„Ž2๐‘ฆโ€ฒโ€ฒ(๐‘ฅ๐‘›) + ๐ถ2,1โ„Ž3๐‘ฆโ€ฒโ€ฒโ€ฒ(๐‘ฅ๐‘›) + โ‹ฏ = 0 ๐ถ0,3

2

๐‘ฆ(๐‘ฅ๐‘›) + ๐ถ13

2,โ„Ž๐‘ฆโ€ฒ(๐‘ฅ๐‘›) + ๐ถ3

2,3

2

โ„Ž2๐‘ฆโ€ฒโ€ฒ(๐‘ฅ๐‘›) + ๐ถ2,3

2

โ„Ž3๐‘ฆโ€ฒโ€ฒโ€ฒ(๐‘ฅ๐‘›) + ๐ถ5

2,3

2

โ„Ž4๐‘ฆโ€ฒ๐‘ฃ(๐‘ฅ๐‘›) + โ‹ฏ = 0 ๐ถ0,3

2

๐‘ฆ(๐‘ฅ๐‘›) + ๐ถ13

2,โ„Ž๐‘ฆโ€ฒ(๐‘ฅ๐‘›) + ๐ถ3 2,3

2

โ„Ž2๐‘ฆโ€ฒโ€ฒ(๐‘ฅ๐‘›) + ๐ถ2,3 2

โ„Ž3๐‘ฆโ€ฒโ€ฒโ€ฒ(๐‘ฅ๐‘›) + ๐ถ5 2,3

2

โ„Ž4๐‘ฆโ€ฒ๐‘ฃ(๐‘ฅ๐‘›) + โ‹ฏ = 0 }

(5)

where (5) is evaluated respectively as follows ๐ถ0,1

2

= ๐›ผ0,1

2

+ ๐›ผ1,1

2

+ ๐›ผ3

2,1

2

= 0 ๐ถ1,1

2

= โˆ’๐›ผ0,1

2

+1

2๐›ผ3

2,1

2

โˆ’ 2๐›ฝ1

2,1

2

= 0 ๐ถ3

2,1

2

= โˆ’1

6๐›ผ0,1

2

+ 1

48๐›ผ3

2,1

2

โˆ’ 1

16๐›ฝ1

2,1

2

= 0}

(6)

๐ถ0,1= ๐›ผ0,1+ ๐›ผ1,1+ ๐›ผ3

2,1+ ๐›ผ2,1= 0 ๐ถ1,1= โˆ’๐›ผ0,1+12๐›ผ3

2,1+ ๐›ผ2,1โˆ’ 2๐›ฝ1,1= 0 ๐ถ3

2,1=1

2๐›ผ0,1+1

8๐›ผ3

2,1+1

2๐›ผ2,1โˆ’ ๐›ฝ1,1= 0 ๐ถ2,1= โˆ’1

6๐›ผ0,1+ 1

48๐›ผ3

2,1+1

6๐›ผ2,1โˆ’1

2๐›ฝ1,1= 0 }

(7)

(5)

๐ถ0,3

2

= ๐›ผ0,3

2

+ ๐›ผ1,3

2

+ ๐›ผ3

2,3

2

+ ๐›ผ2,3

2

+ ๐›ผ5

2,3

2

= 0 ๐ถ1,3

2

= โˆ’๐›ผ0.3

2

+1

2๐›ผ3

2,3

2

+ ๐›ผ2,3

2

+3

2๐›ผ5 2,3

2

โˆ’ 2๐›ฝ3

2,3

2

= 0 ๐ถ3

2,3

2

=1

2๐›ผ0,3

2

+1

8๐›ผ3

2,3

2

+1

2๐›ผ2,3

2

+9

8๐›ผ5

2,3

2

โˆ’ 2๐›ฝ3

2,3

2

= 0 ๐ถ2,3

2

= โˆ’1

6๐›ผ0,3 2

+ 1

48๐›ผ3 2,3

2

+1

6๐›ผ2,3 2

+ 9

16๐›ผ5 2,3

2

โˆ’10

8 ๐›ฝ3 2,3

2

= 0 ๐ถ5

2,3

2

= 1

24๐›ผ0,3

2

+ 1

384๐›ผ3

2,3

2

+ 1

24๐›ผ2,3

2

+ 27

128๐›ผ5

2,3

2

โˆ’ 7

12๐›ฝ3

2,3

2

= 0}

(8)

and

๐ถ0,2= ๐›ผ0,2+ ๐›ผ1,2+ ๐›ผ3

2,2+ ๐›ผ2,2+ ๐›ผ5

2,2+ ๐›ผ3,2= 0 ๐ถ1,2= โˆ’๐›ผ0,2+1

2๐›ผ3

2,2+ ๐›ผ2,2+3

2๐›ผ5

2,2+ 2๐›ผ3,2โˆ’ 2๐›ฝ2,2 = 0 ๐ถ3

2,2=1

2๐›ผ0,2+1

8๐›ผ3 2,2+1

2๐›ผ2,2+9

8๐›ผ5

2,2+ 2๐›ผ3,2โˆ’ 3๐›ฝ2,2= 0 ๐ถ2,2= โˆ’1

6๐›ผ0,2+ 1

48๐›ผ3

2,2+1

6๐›ผ2,2+ 9

16๐›ผ5 2,2+4

3๐›ผ3,2โˆ’5

2๐›ฝ2,2= 0 ๐ถ5

2,2= 1

24๐›ผ0,2+ 1

384๐›ผ3

2,2+ 1

24๐›ผ2,2+ 27

128๐›ผ5

2,2+2

3๐›ผ3,2โˆ’3

2๐›ฝ2,2= 0 ๐ถ3,2= โˆ’1201 ๐›ผ0,2+38401 ๐›ผ3

2,2+1201 ๐›ผ2,2+128081 ๐›ผ5

2,2+154 ๐›ผ3,2โˆ’4572๐›ฝ2,2= 0}

(9)

Normalizing the coefficients ๐›ผ3

2,1

2

, ๐›ผ2,1, ๐›ผ5 2,3

2

& ๐›ผ3,2 of ๐‘ฆ

๐‘›+1

2

, ๐‘ฆ๐‘›+1 , ๐‘ฆ๐‘›+3

2

and ๐‘ฆ๐‘›+2 respectively to 1, solving equation (6), (7), (8) and (9) with the aids of Maple Software for the coefficients of ๐›ผ๐‘—,๐‘– and ๐›ฝ๐‘—,๐‘– and substituting them in (4) gives the first, second, third and fourth point as

๐‘ฆ๐‘›+1

2

= โˆ’1

4๐‘ฆ๐‘›โˆ’1+5

4๐‘ฆ๐‘›+1

8โ„Ž๐‘“๐‘›+1

2

+1

8โ„Ž๐‘“๐‘›โˆ’1

2

๐‘ฆ๐‘›+1=19๐‘ฆ๐‘›โˆ’1โˆ’ 0๐‘ฆ๐‘›+89๐‘ฆ๐‘›+1

2

+13โ„Ž๐‘“๐‘›+1+13โ„Ž๐‘“๐‘›

๐‘ฆ๐‘›+3

2

= โˆ’163

257๐‘ฆ๐‘›โˆ’1+1472

257 ๐‘ฆ๐‘›โˆ’3023

257 ๐‘ฆ๐‘›+1

2

+1971

257 ๐‘ฆ๐‘›+1+273

514โ„Ž๐‘“๐‘›+3

2

+273

514โ„Ž๐‘“๐‘›+1

2

๐‘ฆ๐‘›+2= โˆ’11

65๐‘ฆ๐‘›โˆ’1+29

13๐‘ฆ๐‘›โˆ’64

13๐‘ฆ๐‘›+1

2

+63

13๐‘ฆ๐‘›+1โˆ’64

65๐‘ฆ๐‘›+3

2

+ 6

13โ„Ž๐‘“๐‘›+2+ 6

13โ„Ž๐‘“๐‘›+1 }

(10)

3 ANALYSIS OF THE METHOD

In this section, order and stability properties of the proposed method (10) will be analyzed.

3.1 Order of the Method

In this section, the order of the proposed method (10) will be derived. The method could be transformed to a general matrix form as

(6)

โˆ‘1๐‘—=0๐ถ๐‘—โˆ—๐‘Œ๐‘š+๐‘—โˆ’1= โ„Ž โˆ‘1๐‘—=0๐ท๐‘—โˆ—๐‘Œ๐‘š+๐‘—โˆ’1, (11) where C and D are constant coefficient matrices of the method (11) is equivalent to the following form

โˆ’1

4๐‘ฆ๐‘›โˆ’1+ 5

4๐‘ฆ๐‘› + ๐‘ฆ๐‘›+1

2

= โˆ’1

8โ„Ž๐‘“๐‘›+1

2

โˆ’1

8โ„Ž๐‘“๐‘›โˆ’1

2

1

9๐‘ฆ๐‘›โˆ’1โˆ’ 0๐‘ฆ๐‘›+8 9๐‘ฆ

๐‘›+1 2

+ ๐‘ฆ๐‘›+1 = โˆ’1

3โ„Ž๐‘“๐‘›+1โˆ’1 3โ„Ž๐‘“๐‘›

โˆ’163

257๐‘ฆ๐‘›โˆ’1+1472

257 ๐‘ฆ๐‘›โˆ’3023

257 ๐‘ฆ๐‘›+1

2

+1971

257 ๐‘ฆ๐‘›+1+ ๐‘ฆ๐‘›+3

2

= โˆ’273

514โ„Ž๐‘“๐‘›+3

2

โˆ’273

514โ„Ž๐‘“๐‘›+1

2

(12)

โˆ’11

65๐‘ฆ๐‘›โˆ’1+29

13๐‘ฆ๐‘›โˆ’64 13๐‘ฆ

๐‘›+1 2

+63

13๐‘ฆ๐‘›+1โˆ’64 65๐‘ฆ

๐‘›+3 2

+ ๐‘ฆ๐‘›+2= โˆ’ 6

13โ„Ž๐‘“๐‘›+2โˆ’ 6 13โ„Ž๐‘“๐‘›+1

also (12) can be written as

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

โˆ’

โˆ’

โˆ’

13 0 29 65 0 11

257 0 1472 257

0 163

0 9 0

0 1

4 0 5 4 0 1

[ ๐‘ฆ๐‘›โˆ’3 ๐‘ฆ๐‘›โˆ’12

๐‘ฆ๐‘›โˆ’1 ๐‘ฆ๐‘›2]

+

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช

๏ƒซ

๏ƒฉ

โˆ’

โˆ’

โˆ’

65 1 64 13

63 13

64

0 257 1

1971 257

3023

0 0 9 1

8

0 0 0

1

[ ๐‘ฆ๐‘›+1 ๐‘ฆ๐‘›+12

๐‘ฆ๐‘›+3 ๐‘ฆ๐‘›+22]

3 2

1 2

1

0 0 1 0

8 0 0 0 1

3

0 0 0 0

0 0 0 0

n

n n n

f h f

f f

โˆ’

โˆ’

โˆ’

๏ƒฉ โˆ’ ๏ƒน๏ƒฉ ๏ƒน

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช โˆ’ ๏ƒบ

= ๏ƒช ๏ƒบ๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ ๏ƒซ ๏ƒป

๏ƒช ๏ƒบ

๏ƒซ ๏ƒป

+h

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

โˆ’

โˆ’

โˆ’

โˆ’

โˆ’

โˆ’

13 0 6

13 0 6

514 0 0 273

514 273

0 3 0

0 1

0 0

8 0 1

[ ๐‘“๐‘›+1

2

๐‘“๐‘›+1

๐‘“๐‘›+3

2

๐‘“๐‘›+2]

(13)

where

(7)

๐ถ0=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

0 0 0 0

๐ถ1=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

โˆ’

โˆ’

โˆ’

65 11 257 163

9 14

1

๐ถ2=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

0 0 0 0

๐ถ3=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช

๏ƒซ

๏ƒฉ

13 29 257 1472

0 4 5

๐ถ4=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช

๏ƒซ

๏ƒฉ

โˆ’

13 64 257 3023

9 8 1

๐ถ5=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

13 63 257 1971

1 0

๐ถ6 =

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช

๏ƒซ

๏ƒฉ

โˆ’13 64 1 0 0

๐ถ7=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

1 0 0 0

๐ท0=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

0 0 0 0

๐ท1=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

0 0 0 0

๐ท2=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช

๏ƒซ

๏ƒฉโˆ’

0 0 0 8 1

๐ท3=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช

๏ƒซ

๏ƒฉ

โˆ’

0 0 3 1 0

๐ท4=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

โˆ’

โˆ’

0 514 273 08

1

๐ท5=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช

๏ƒซ

๏ƒฉ

โˆ’

โˆ’

13 6 0

3 1 0

๐ท6=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช

๏ƒซ

๏ƒฉ

โˆ’ 0 514 273

0 0

๐ท7=

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช

๏ƒซ

๏ƒฉ

โˆ’13 6 0 0 0

Definition 3.1.1 According to [26], the order of the block method (10) and its associated linear operator are given by

7 7

0 0

[ ( ); ] [ j ( )] j '( )]

j j

L y x h c y x jh h D y x jh

= =

=

๏ƒฅ

+ โˆ’

๏ƒฅ

+ (14)

where p is unique integer such that Eq = 0, q = 0,1, โ€ฆ p and Ep+1โ‰  0, where the Eq are constant matrix with

๐ธ0 = โˆ‘ ๐ถ๐ฝ

7 ๐ฝ=0

= 0

๐ธ1= โˆ‘[๐‘—๐ถ๐‘—โˆ’ 2๐ท๐‘—]

7

๐ฝ=0

= 0

๐ธ2= โˆ‘ [1

2!๐‘—2๐ถ๐‘—โˆ’ 2๐‘—๐ท๐‘—]

7

๐ฝ=0

= 0

๐ธ3 = โˆ‘ [1

3!๐‘—3๐ถ๐‘—โˆ’ 21 2!๐‘—2๐ท๐‘—]

7

๐ฝ=0

= 0

๐ธ4 = โˆ‘ [1

4!๐‘—4๐ถ๐‘—โˆ’ 21 3!๐‘—3๐ท๐‘—]

7

๐ฝ=0

= 0

(8)

๐ธ5 = โˆ‘ [1

5!๐‘—5๐ถ๐‘—โˆ’ 21 4!๐‘—4๐ท๐‘—]

7

๐ฝ=0

= 0

๐ธ6 = โˆ‘ [1

6!๐‘—6๐ถ๐‘—โˆ’ 21 5!๐‘—5๐ท๐‘—]

7

๐ฝ=0

= 0

๐ธ7 = โˆ‘ [1

7!๐‘—7๐ถ๐‘—โˆ’ 21

6!๐‘—6๐ท๐‘—]

7๐ฝ=0 =

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง

๏ƒจ

๏ƒฆ

โˆ’

โˆ’

โˆ’

โˆ’

5123 192 7385

274 5065

221 9353

121

๏‚น

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

0 0 0 0

Therefore, the developed method is of order 6, with error constant

E7=

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง

๏ƒจ

๏ƒฆ

โˆ’

โˆ’

โˆ’

โˆ’

5123 192 7385

274 5065

221 9353

121

(15)

3.2 Stability Analysis of the Method

In this section, we investigate the Zero and A- stability property of the proposed method (10).

Definition 3.2.1 A linear multistep method is said to be zero stable if no root of the first characteristics polynomial has modulus higher than 1 and that any root with modulus 1 is simple [31].

Definition 3.2.2 A linear multistep method is said to be an A-stable method if its region of stability encloses the entire negative half-plane [31].

The stability of the scheme (10-11) can be obtained by applying the standard test equation of the form

๐‘ฆโ€ฒ = สŽ๐‘ฆ, ๐‘…๐‘’(สŽ) < 0 (16) where สŽ is a complex number.

To get the following solutions

(9)

๐‘ฆ๐‘›+1

2

= โˆ’1

4๐‘ฆ๐‘›โˆ’1+5

4๐‘ฆ๐‘›+1

8โ„ŽสŽ๐‘ฆ๐‘›+1

2

+1

8โ„ŽสŽ๐‘ฆ๐‘›โˆ’1

2

๐‘ฆ๐‘›+1=1

9๐‘ฆ๐‘›โˆ’1โˆ’ 0๐‘ฆ๐‘›+8

9๐‘ฆ๐‘›+1 2

+1

3โ„ŽสŽ๐‘ฆ๐‘›+1+1

3โ„ŽสŽ๐‘ฆ๐‘› ๐‘ฆ๐‘›+3

2

= โˆ’163257๐‘ฆ๐‘›โˆ’1+1472257 ๐‘ฆ๐‘›โˆ’3023257 ๐‘ฆ๐‘›+1

2

+1971257 ๐‘ฆ๐‘›+1+273514โ„ŽสŽ๐‘ฆ๐‘›+3

2

+273514โ„ŽสŽ๐‘ฆ๐‘›+1

2

๐‘ฆ๐‘›+2= โˆ’11

65๐‘ฆ๐‘›โˆ’1+29

13๐‘ฆ๐‘›โˆ’64

13๐‘ฆ๐‘›+1

2

+63

13๐‘ฆ๐‘›+1โˆ’64

65๐‘ฆ๐‘›+3

2

+ 6

13โ„ŽสŽ๐‘ฆ๐‘›+2+ 6

13โ„ŽสŽ๐‘ฆ๐‘›+1 }

(17)

(17) can also be written as

1 1 0 0 0

8

8 1

1 0 0

9 3

3023 273 1971 273

1 0

257 514 257 514

64 63 6 64 6

13 13 13 65 1 13

h

h

h h

h h

๏ฌ

๏ฌ

๏ฌ ๏ฌ

๏ฌ ๏ฌ

๏ƒฉ โˆ’ ๏ƒน

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช โˆ’ โˆ’ ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช โˆ’ โˆ’ โˆ’ ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช โˆ’ โˆ’ โˆ’ ๏ƒบ

๏ƒซ ๏ƒป

[ ๐‘ฆ๐‘›+1 ๐‘ฆ๐‘›+12

๐‘ฆ๐‘›+3 ๐‘ฆ๐‘›+22]

=

1 1 5

0 4 8 4

1 1

0 0

9 3

163 1472

0 0

257 257

11 29

0 0

65 13

h h

๏ฌ

๏ฌ

๏ƒฉ โˆ’ ๏ƒน

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช โˆ’ ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช ๏ƒบ

๏ƒช โˆ’ ๏ƒบ

๏ƒซ ๏ƒป

[ ๐‘ฆ๐‘›โˆ’3 ๐‘ฆ๐‘›โˆ’12

๐‘ฆ๐‘›โˆ’1 ๐‘ฆ๐‘›2]

(18)

From (18) it is given that

๐ด๐‘Œ๐‘š = ๐ต๐‘Œ๐‘šโˆ’1 (19)

If m is the number of block and r is the number of points in the block, then n = mr Here, r = 2 and n = 2m. It follows that

๐‘Œ๐‘š = [

๐‘ฆ2๐‘š+1

๐‘ฆ2๐‘š+12

๐‘ฆ2๐‘š++3

๐‘ฆ3๐‘š+22]

= [

๐‘ฆ๐‘›+1

๐‘ฆ๐‘›+12

๐‘ฆ๐‘›+3

๐‘ฆ๐‘›+222]

, ๐‘Œ๐‘šโˆ’1= [

๐‘ฆ2(๐‘šโˆ’1)โˆ’3

๐‘ฆ2(๐‘šโˆ’1)โˆ’12

๐‘ฆ2(๐‘šโˆ’1)โˆ’1

๐‘ฆ2(๐‘šโˆ’1)2]

= [

๐‘ฆ๐‘›โˆ’3

๐‘ฆ๐‘›โˆ’12

๐‘ฆ๐‘›โˆ’1

๐‘ฆ๐‘›2]

and the coefficient matrices are given as

๐ด =

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

โˆ’

โˆ’

โˆ’

โˆ’

โˆ’

โˆ’

โˆ’

โˆ’

โˆ’

๏ฌ

๏ฌ

๏ฌ

๏ฌ

๏ฌ

๏ฌ

h h

h h

h h

13 1 6 65

64 13

6 13 63 13

64

514 0 1 273 257

1971 514

273 257

3023

0 3 0

1 1 9

8

0 0

8 0 1 1

(10)

๐ต =

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒบ๏ƒบ

๏ƒป

๏ƒน

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒช๏ƒช

๏ƒซ

๏ƒฉ

โˆ’

โˆ’

โˆ’

13 0 29 65

0 11

257 0 1472 257

0 163

3 0 1 9

0 1

4 5 8

1 4 0 1

๏ฌ

๏ฌ h h

The stability polynomial of the proposed method was computed with the aid of Maple Software and the result is found to be

๐‘‘๐‘’๐‘ก(๐ด๐‘ก โˆ’ ๐ต) =

(20)

๐‘…(๐‘ก, 0) = (21)

๐‘ก = 0 , 0 , 1, โˆ’0.2629880097

3.3 A- stability of the Proposed Method

In this section, the region for the absolute stability of the proposed methods is plotted, by considering the stability polynomials (20). The set of point defined by ๐‘ก = ๐‘’๐‘–๐›ณ, 0 โ‰ค ๐œƒ โ‰ค 2๐œ‹ describes the boundary of the stability region. The following stability region was the complex plot of the proposed method with the aid of Maple Software.

(11)

Figure 2: A-stability region of the Proposed Method (RDIBM)

4 IMPLEMENTATION OF THE METHOD

Consider the system of first order initial value problem of ordinary differential equation of the form

๐‘ฆโ€ฒ = ๐‘“(๐‘ฅ, ๐‘Œฬ‚), ๐‘Œฬ‚(๐‘Ž) = ีฒ, ๐‘Ž โ‰ค ๐‘ฅ โ‰ค ๐‘ (22)

๐‘Œฬ‚ = (๐‘ฆ1, ๐‘ฆ2, ๐‘ฆ3, โ€ฆ โ€ฆ โ€ฆ ๐‘ฆ๐‘›), (23)

using Newtonโ€™s iteration to implement the methods (10).

Let ๐‘ฆ๐‘– and ๐‘ฆ(๐‘ฅ๐‘–) be the approximate and exact solutions of system (22-23) respectively.

Define the absolute error as

(๐‘’๐‘Ÿ๐‘Ÿ๐‘œ๐‘Ÿ๐‘–)๐‘ก = |(๐‘ฆ๐‘–)๐‘กโˆ’ (๐‘ฆ(๐‘ฅ๐‘–))๐‘ก| (24)

and the maximum error as MAXE = ๐‘š๐‘Ž๐‘ฅโŸ

1โ‰ค๐‘–โ‰ค๐œ

(๐‘š๐‘Ž๐‘ฅ (๐‘’๐‘Ÿ๐‘Ÿ๐‘œ๐‘ŸโŸ ๐‘–)๐‘ก

1โ‰ค๐‘–โ‰คโ„ต

) (25)

๐œ and โ„ต are the total number of step and equations respectively.

From the method (10)

Im

Re stable

stable stable

stable

unstable

unstable

(12)

1 1 1 1 1

2 2 2

2 1 1 1 2

2

3 3 1 1 3 1 3

2 2 2 2

4 2 1 1 3 2 1 4

2 2

1 1

8 8

8 1 1

9 3 3

3023 1971 273 273

257 257 514 514

64 63 64 6 6

13 13 65 13 13

n n n

n n n

n

n n n n n

n n n n

n n

F y hf hf

F y y hf hf

F y y y hf hf

F y y y y hf hf

๏ฅ

๏ฅ

๏ฅ

๏ฅ

+ + โˆ’

+ + +

+ + + + +

+ + + + + +

= โˆ’ โˆ’ โˆ’

= โˆ’ โˆ’ โˆ’ โˆ’

= + โˆ’ โˆ’ โˆ’ โˆ’

= + โˆ’ + โˆ’ โˆ’ โˆ’

(26)

The โ„ฐ1, โ„ฐ2, โ„ฐ2 and โ„ฐ4 are the back values of (10) as

1 1

2 1

3 1

4 1

1 5

4 4

1 0

9

163 1472

257 257

11 29

65 13

n n

n n

n n

n n

y y

y y

y y

y y

๏ฅ

๏ฅ

๏ฅ

๏ฅ

โˆ’

โˆ’

โˆ’

โˆ’

= โˆ’ +

= โˆ’

= โˆ’ +

= โˆ’ +

(27)

Let ๐‘ฆ๐‘›+๐‘—(๐‘–+1), ๐‘— =1

2, 1,3

2, 2, denote the (๐‘– + 1)๐‘กโ„Ž iterative values of ๐‘ฆ๐‘›+๐‘— and consider ๐‘’๐‘›+๐‘—(๐‘–+1) = ๐‘ฆ๐‘›+๐‘—(๐‘–+1)โˆ’ ๐‘ฆ๐‘›+๐‘—(๐‘–) , ๐‘— =1

2, 1,3

2, 2 (28)

Now, the Newtonโ€™s iteration for the proposed method will have the form ๐‘ฆ๐‘›+๐‘—(๐‘–+1) = ๐‘ฆ๐‘›+๐‘—(๐‘–) โˆ’ (๐น๐‘—(๐‘ฆ๐‘›+๐‘—

(๐‘–) )) (๐นโ€ฒ๐‘—(๐‘ฆ๐‘›+๐‘—(๐‘–) )), ๐‘— =1

2, 1,3

2, 2 (29)

๐‘ฆ๐‘›+๐‘—(๐‘–+1) = ๐‘ฆ๐‘›+๐‘—(๐‘–) โˆ’ (๐นโ€ฒ๐‘—(๐‘ฆ๐‘›+๐‘—(๐‘–) ))

โˆ’1

(๐น๐‘—(๐‘ฆ๐‘›+๐‘—(๐‘–) )) , ๐‘— =1

2, 1,3

2, 2 (30)

๐‘ฆ๐‘›+๐‘—(๐‘–+1)โˆ’ ๐‘ฆ๐‘›+๐‘—(๐‘–) = โˆ’ (๐นโ€ฒ๐‘—(๐‘ฆ๐‘›+๐‘—(๐‘–) ))

โˆ’1

(๐น๐‘—(๐‘ฆ๐‘›+๐‘—(๐‘–) )) , ๐‘— =1

2, 1,3

2, 2 (31)

๐‘’๐‘›+๐‘—(๐‘–+1) = โˆ’ (๐นโ€ฒ๐‘—(๐‘ฆ๐‘›+๐‘—(๐‘–) ))

โˆ’1

(๐น๐‘—(๐‘ฆ๐‘›+๐‘—(๐‘–) )) , ๐‘— =1

2, 1,3

2, 2 (32)

It can be written as

(13)

(๐นโ€ฒ๐‘—(๐‘ฆ๐‘›+๐‘—(๐‘–) )) ๐‘’๐‘›+๐‘—(๐‘–+1)= โˆ’ (๐น๐‘—(๐‘ฆ๐‘›+๐‘—(๐‘–) )) , ๐‘— =1

2, 1,3

2, 2 (33)

Equation (33) will also be written in its matrix form as:

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด ๏€ณ

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด

๏€ด ๏€ฒ

๏€ฑ

trix jacobianm a

i n

i n i

n i n

i n

i n i

n i n

i n

i n i

n i n

y h F y

F

y F y

F

y F y

F

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒจ

๏ƒฆ

๏‚ถ

โˆ’ ๏‚ถ

๏‚ถ

โˆ’ ๏‚ถ

โˆ’

๏‚ถ

๏‚ถ

โˆ’

๏‚ถ โˆ’

๏‚ถ

โˆ’

๏‚ถ

โˆ’ ๏‚ถ

โˆ’

๏‚ถ

๏‚ถ

โˆ’

+ + +

+

+ +

+ +

+ + +

+

) (

2 ) (

2 )

( 1 ) (

1

) (

2 3 ) (

2 3 )

( 2 1 ) (

2 1

) (

1 ) (

1 )

( 2 1 ) (

2 1

13 1 6 65

64 13

6 13 63 13

64

514 0 1 273 257

1971 514

273 257 3023

0 3 0

1 1 9

8

0 0

8 0 1 1

=

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒจ

๏ƒฆ

++ + + ++ + +

) 1 (

2 ) 1 (

2 3 ) 1 (

1 ) 1 (

2 1

i n

i n

i n

i n

e e e e

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง

๏ƒจ

๏ƒฆ

โˆ’

โˆ’

โˆ’

โˆ’

โˆ’

โˆ’

โˆ’

65 1 64 13

63 13

64

0 257 1

1971 257

3023

0 0 9 1

81 0 0 0

+

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒจ

๏ƒฆ

0 0 0 0

0 0 0

0 3

0 1 0 0

8 0 0 1 0

h +

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง

๏ƒจ

๏ƒฆ

13 0 6 13 0 6

514 0 0 273 514 273

0 3 0

0 1

0 0 8 0

1

h

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒจ

๏ƒฆ

+ + + +

2 2 3 1 2 1

n n n n

f f f f

+

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒจ

๏ƒฆ

2 2 3 1 2 1

๏ฅ

๏ฅ

๏ฅ

๏ฅ

(34)

4.1 Test Problems

To validate the method developed, (RDIBM), a code in โ€˜Cโ€™ (programming Language) with Equation (34) would be used to solve the following stiff IVPs.

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒจ

๏ƒฆ

+ + + +

2 2 3 1 2 1

n n n n

y y y y

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒท๏ƒท

๏ƒธ

๏ƒถ

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒง๏ƒง

๏ƒจ

๏ƒฆ

โˆ’

โˆ’

โˆ’

n n n n

f f

f f

2 1 1 2 3

(14)

Table 1: Sample of First Order Initial Value Problem of Stiff ODEs S/

n

Problems Initial Condition s

Interval Exact Solutions Eigen Values

Sou rce 1 ๐‘ฆ1โ€ฒ = 198y1+ 199y1

๐‘ฆ2โ€ฒ

= โˆ’398๐‘ฆ1โˆ’ 399๐‘ฆ2

y1(0) = 1 y2(0)

= โˆ’1

0 โ‰ค ๐‘ฅ โ‰ค

10 ๐‘ฆ1(๐‘ฅ) = ๐‘’โˆ’๐‘ฅ ๐‘ฆ2(๐‘ฅ) = โˆ’๐‘’โˆ’๐‘ฅ

โˆ’1,โˆ’20

0 [4]

2 ๐‘ฆ1โ€ฒ = y2 ๐‘ฆ2โ€ฒ = โˆ’2๐‘ฆ1 ๐‘ฆ3โ€ฒ = y2+ 2๐‘ฆ3

y1(0) = 0 y2(0) = 0 y3(0) = 1

0 โ‰ค ๐‘ฅ โ‰ค

4๐œ‹ ๐‘ฆ1(๐‘ฅ) = 2๐‘๐‘œ๐‘ ๐‘ฅ + 6๐‘ ๐‘–๐‘›๐‘ฅ

โˆ’ 6๐‘ฅ โˆ’ 2 ๐‘ฆ2(๐‘ฅ) = โˆ’2๐‘ ๐‘–๐‘›๐‘ฅ + 6๐‘๐‘œ๐‘ ๐‘ฅ

โˆ’ 6 ๐‘ฆ3(๐‘ฅ) = 2๐‘ ๐‘–๐‘›๐‘ฅ โˆ’ 2๐‘๐‘œ๐‘ ๐‘ฅ

+ 3

[25]

3 ๐‘ฆโ€ฒ = 5๐‘’5๐‘ฅ(๐‘ฆ โˆ’ ๐‘ฅ)2 + 1

y(0) = 0 0 โ‰ค ๐‘ฅ โ‰ค 1 ๐‘ฆ1(๐‘ฅ) = ๐‘ฅ โˆ’ ๐‘’โˆ’5๐‘ฅ [4]

4 ๐‘ฆ1โ€ฒ = โˆ’20๐‘ฆ1โˆ’19y2

y2โ€ฒ = โˆ’19y1โˆ’20y2 y1(0) = 2 y2(0) = 0

0 โ‰ค ๐‘ฅ โ‰ค

20 ๐‘ฆ1(๐‘ฅ) = ๐‘’โˆ’๐‘ฅ ๐‘ฆ2(๐‘ฅ) = โˆ’๐‘’โˆ’๐‘ฅ

โˆ’1,โˆ’20 0

[5]

5 RESULT AND DISCUSSIONS

The sample problems presented in this paper are solved using the proposed methods. The result of the tested problems are tabulated and compared with the existing ones. The graphs highlighting the performance of these methods are plotted. The acronyms below are used in the tables.

h= step-size;

MHTD =Method

MAX-ERR = Maximum Error;

EXEC-TIME= Executional Time in second;

ABISBDF = An A-stable Block Integrator Scheme for the Solution of First Order System of IVP of Ordinary Differential Equations

3ESBBDF = Extended 3-Point Super class of Block Backward Differentiation formula for Solving Stiff Initial Value Problems.

RDIBM = A Robust Diagonally Implicit Block Method with Two Off-Step Points for Solving First Order Stiff IVP of ODEs

3NBBDF = A New Fifth Order implicit block method for Solving First Order Stiff Ordinary Differential Equations

3BDF = Implicit r-point block backward differentiation formula for solving first-order stiff ODEs

(15)

Table 2: Comparison of Errors for Problem 1

Figure 3: Graph of Log10(๐‘€๐ด๐‘‹๐ธ) against the step size h for Problem 1

0 0.001 0.002 0.003 0.004 0.005 0.006 0.007 0.008 0.009 0.01

-9 -8 -7 -6 -5 -4 -3 -2 -1

h

logError

3BBDF ABISBDF RDOBBDF 3NBBDF

Numerical Result for Problem 1

๐’‰ MTHD NS MAX-ERR EXEC-TIME 10โˆ’2 3NBBDF

ABISBDF RDIBM 3BDF

333 555 100 333

1.94447e-004 5.83217e-003 1.52564e-004 1.07308e-02

1.20394e-002 5.68676e-002 3.93719e-003 31,867ฮผs 10โˆ’3 3NBBDF

ABISBDF RDIBM 3BDF

3,333 5,555 1,000 3,333

2.07993e-006 6.05338e-005 1.76763e-006 1.10060e-03

1.19193e-001 5.64515e-001 1.87573e-002 258,361ฮผs 10โˆ’4 3NBBDF

ABISBDF RDIBM 3BDF

33,333 55,555 10,000 33,333

2.09995e-008 6.26692e-007 1.79766e-008 1.10333e-04

1.19296e+000 5.68143e+000 1.66571e-001 2,582,756ฮผs 10โˆ’5 3NBBDF

ABISBDF RDIBM 3BDF

333,333 555,555 100,000 333,333

2.10257e-010 6.32740e-009 1.82566e-010 1.10361e-05

1.19173e+001 5.59821e+001 1.43458e+000 26,011,417ฮผs 10โˆ’6 3NBBDF

ABISBDF RDIBM 3BDF

3,333,333 5,555555 1,000,000 3,333,333

1.41029e-011 6.33362e-011 1.85567e-012 1.10363e-06

1.19110e+002 5.53567e+002 1.28786e+001 260,435,329ฮผs

(16)

Table 3 Comparison of Errors for Problem 2

Figure 4: Graph of Log10(๐‘€๐ด๐‘‹๐ธ) against the step size h for Problem 2

0 0.001 0.002 0.003 0.004 0.005 0.006 0.007 0.008 0.009 0.01

-9 -8 -7 -6 -5 -4 -3 -2 -1

h

logError

3ESBBDF ABISBDF RDOBBDF

Numerical Result for Problem 2

๐‘ฏ MTHD NS MAX-ERR EXEC-TIME

10โˆ’2 3ESBBDF ABISBDF RDIBM

333 555 100

2.5780e-002 3.8321e-003 2.2979e-003

3.96563e-001 5.58677e-002 2.98278e-002 10โˆ’3 3ESBBDF

ABISBDF RDIBM

3,333 5,555 1,000

2.2907e-003 4.0533e-005 2.1066e-005

3.66799e+000 5.54512e-001 1.79012e-001 10โˆ’4 3ESBBDF

ABISBDF RDIBM

33,333 55,555 10,000

2.0972e-005 4.2669e-007 1.9765e-007

3.35906e+000 5.52149e-001 1.49191e+000 10โˆ’5 3ESBBDF

ABISBDF RDIBM

333,333 555,555 100,000

2.002e-007 4.3274e-009 1.7932e-009

3.01024e+001 5.49867e+000 1.20674e+001 10โˆ’6 3ESBBDF

ABISBDF RDIBM

3,333,333 5,555555 1,000,000

1.8702e-009 4.3335e-009 1.6004e-011

2.96155e+002 5.35204e+000 1.20678e+001

Figure

Updating...

References

Related subjects :