297 American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) ISSN (Print) 2313-4410, ISSN (Online) 2313-4402 © Global Society of Scientific Research and Researchers http://asrjetsjournal.org/ Implicit Second Derivative Hybrid Linear Multistep Method with Nested Predictors for Ordinary Differential Equations S. E. Ekoroa*, M. N. O. Ikhileb, I. M. Esuabanac aDepartment of Mathematics, University of Calabar, Calabar, Cross River State, Nigeria bDepartment of Mathematics, University of Benin, Benin, Edo State, Nigeria cDepartment of Mathematics, University of Calabar, Calabar, Cross River State, Nigeria aEmail: ekorosam@yahoo.com bEmail: mnoikhile@yahoo.com cEmail: esuabanaita@gmail.com Abstract In this paper, we considered an implicit hybrid linear multistep method with nested hybrid predictors for solving first order initial value problems in ordinary differential equations. The derivation of the methods is based on interpolation and collocation approach using polynomial basis function. The region of absolute stability of the method is investigated using the boundary locus approach and the methods have been found to be A − stable for step-length 6.k ≤ Keywords: Linear multistep methods; hybrid; nesting; interpolation; collocation; boundary locus. 1. Introduction The conventional linear multistep method (LMM) is defined as 0 0 k k j n j j n j j j y h fα β+ + = = =∑ ∑ (1.1) ------------------------------------------------------------------------ * Corresponding author. http://asrjetsjournal.org/ American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 298 where jα and jβ are parameter constants to be determined. The kβ determines if the linear multistep method is explicit or implicit. For explicit LMM (1.1), 0kβ = and for implicit methods, 0kβ ≠ . This is a popular method for the numerical approximation of the solutions of initial value problems in ordinary differential equations ( )' ,y f x y= , ( )0 0y x y= (1.2) Its stability and order are subject to some constraints by [4]. Modification have been made to overcome the barrier, see [2,5,6,7,15,16] among others. Reference [6] introduced a second derivative term into the Adams- type LMM (1.1) to obtain the second derivative linear multistep (SDLMM) of the form 2 1 1 0 ' k n k k n k j n j n k j y y h f h fα β+ − + − + + = = + +∑ (1.3) Off-step points have been introduced into this linear multistep method to overcome Dahlquist order and stability barrier. Other extension of (1.1) can be found in [10,1,8,3,11,14,16]. Our interest in this paper is to construct an implicit second derivative hybrid linear multistep method of the form ( ) ( ) ( )2 1 0 ' mm k m m m n k n k j n j n v n kkv j y y h f f h fβ β λ+ + − + + + =   = + + +     ∑ (1.4) which are of order 3p k= + with the hybrids ( ) ( ) ( ) 1 2 0 ' l l ll l k l l l n v n k j n j n v n vv v j y y h f f h fβ β λ ++ + + + + =   = + + +     ∑ (1.5) of order * 4p k= + , where ( ) ( ) ( ) 0 2 0 ' k l l l n v j n j n k n kk k j y y h f h fα β λ− − − + + + + = = + +∑ (1.6) of order ** 2p k= + for 0(1)m 1l = − This method (1.4) seeks to approximate the solution of (1.2). The idea is to approximate (1.2) through the integration interval [ ]0 , xNx where ( ) :y x [ ]0 , xNx mℜ→ in which [ ]0: , x m Nf x × ℜ is smooth. 2. Specification of the hybrid methods (1.4) The hybrid methods (1.4) with the hybrid predictors (1.5) and (1.6) have constant parameters American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 299 ( ){ } 0 km j j β = , ( ) m m vβ , ( )m kλ , ( ){ } 0 kl j j β = , ( ) l l vβ , ( ) l l vλ , ( ){ } 0 kl j j λ − = , ( )l kβ − and ( )l kλ − to be determined in such a way that the hybrid method (1.4) become stable. The method (1.4) is the hybrid method of Adams-type equipped with nested functions evaluation of the hybrid predictors (1.5) and (1.6). The hybrid parameters are chosen according as 1 2mv k= − , 1 2 l l v kv + + = , 0(1)m 1,l = − ( )0,lv k∈ , ,lv j≠ 0(1)kj = , 1,2,3,...k = , 1m k= − 2.1 Construction of the Hybrid methods (1.4) We assume the solution of (1.4) of the form 3 0 (x) k j j j y a x + = = ∑ (2.1) where{ } 3 0 k j j a + = are real constant parameters to be determined and ( ), (1)k 3jx j o = + is the polynomial basis function. Differentiating (2.1) twice to obtain ( ) 3 1 1 '(x) , k j j j y f x y ja x + − = = = ∑ (2.2) 3 2 2 ''(x) '(x, y) (j 1) k j j j y f j a x + − = = = −∑ (2.3) Interpolating (2.1), (2.2) and (2.3) at n kx x += and collocating (2.2) at n jx x += , 0(1)k 2j = − and mn vx x += we obtain the system of equations ( ) ( ) ( ) ( ) ( )( ) 2 3 1 1 1 3 32 11 3 3 3 1 . . . 0 1 2 . . . 3 320 1 . . . ... . . ... . . ... . . 20 1 . . . 3 0 1 . . .2 3 0 1 2 . . . 3 2 m m k n k n k n k k n n k nn k n k n k n v n v k n k x x x x k x k xx x k x x k x K k x + + − + − + − + + ++ + + + + + + +     +    +             +   +    + +  0 1 2 1 2 3 . . . k k k a a a a a a + + +                             = 1 1 . . . ' m n k n n n k n v n k y f f f f f + − + + + +                              (2.4) Solving equation (2.4) with MATHEMATICA 10.0 Software package, the coefficients American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 300 ( )( )' 0 1 3ja s j k= + are obtained. Substituting these coefficients into (2.1) yields the discrete scheme for each k . 3. Construction of the hybrid Predictors The corresponding hybrid predictor is obtained from the polynomial interpolant ( ) 4 1 0 k j n l j j y x v h b x + + = + = ∑ (3.1) where { } 4k j j o b + = are parameter constants to be determined, { } 4 0 kj j x + = is the polynomial basis function. Following the approach as in section (3), we obtain the system of equations ( ) ( ) ( ) 2 3 4 2 2 1 1 2 1 1 . . . 0 1 2 . . . 3 0 1 2 . . . 3 . . . . . . . . . . . . . . . 0 0 2 . . . 3 0 0 2 . . . 20 2 k n k n n k k n n k n n k n v n v k n v x x x x k x k x x k x k x + + + + + + + + + + + + +      +   +             +   +  0 1 2 2 3 . . . k k a a a a a + +                          = 1 . . . ' n k n n n v n v y f f f f + + + +                          (3.2) Equation (3.2) is solved with MATHEMATICA 10.0 software package to obtain the coefficients of the hybrid predictor (1.5) The corresponding error constants for the hybrid scheme and its hybrid predictors are obtained for each value of k from the Taylor series expansion of (1.4), (1.5) and (1.6) about nx . These are respectively ( ) ( ) ( )1 1 2 1 0p p p n k n k p ny y x C h y x h+ + + + + +− = + (3.3) ( ) ( ) ( )* * * *1 1 1 1 2 1 0 l l p p p n v n v npy y x C h y x h + + + + + + + +− = + (3.4) ( ) ( )** ** **0 0 1 2 1 (x ) 0p p n v n v npy y x C h h+ + + + +− = + (3.5) where ( )n ky x + , ( )1ln vy x ++ and ( )0n vy x + are the theoretical solutions; 1pC + , * 1pC + and ** 1pC + are error American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 301 constants of (1.4) ,(1.5) and (1.6) respectively. Due to the processing speed and the memory capacity of the laptop computer used in the derivation, only few stable members of the family of the method could be obtained. If the method can be derived using higher processor, more stable members can be obtained from step- number 10k ≥ . Examples of A − stable members of the family of the hybrid methods (1.4) with error constants are: For 0 11, 0,v 2 k m= = = 1 1 1 2 2 6 3 6 n n n nn f fy h f y+ + +    = + + +     , 5 1 2880 C = − with hybrid 21 1 1 1 2 73 1 ' 8 8 8 16 n n n nn y yy hf h f+ + + + = − + + + , 6 1 384 C − = For 1 32, 1,v 2 k m= = = 21 1 2 3 1 2 2 11 4728 1 ' 720 60 45 240 120 n n n n n nn f f fy h f y h f+ + + + + +    = − + + + + −     , 6 1 14400 C = − with hybrids 1 0 3 7,v 2 4 v = = 1 2 3 7 2 7 2 4 4 92516 29 3920 135 6615 80 1260 n n n nn n n f f fy h f y f+ + + + + +    = − − − + +     , 7 7 40960 C = 2 2 1 2 2 7 4 231 3 7 1995 21 ' 1024 2048 256 2048 90 n n n n n n f y y y h fy h + + + + + = − − + + + , 6 1 50812 C = − For 2 53, 2,v , 6 2 K m p= = = = 21 2 3 3 2 5 3 2 23 223136 1 ' 5400 360 120 225 1080 90 n n n n n n nn f f f fy y h f h f+ + + + + + +    = + − + + + −     , 7 13 604800 C = − American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 302 with hybrids 2 1 0 5 11 23, v ,v 2 4 8 v = = = 21 2 5 11 3 3 11 2 4 4 13 3 167 44048 203 2 ' 232320 4480 17280 144345 1920 99 n n n n nn n n f f fy h f f y h f+ + + + + + +    = − + − − − + +     8 1 107520 C = 2 1 2 1 11 23 3 23 4 8 8 77 359 269 40716032868 286 ' 16250880 6912000 501760 87483375 6144 36225 n n n n nn n n f f f f hy h f y f+ + + + + + +    = − + − − − + +     8 1513 1651507200 C = 2 3 1 2 3 3 23 8 47495 35 161 345 2348185 805 ' 393216 589824 262144 65536 2359296 131072 n n n n n n n hf y y y y h fy + + + + + + = − + − + + + , 7 161 12582912 C = 4. Stability of the Hybrid Schemes (1.4) This section considers some important definitions and stability properties of the hybrid schemes. Definition 1: A numerical scheme (1.4) is A − stable if the region of absolute stability lies entirely in the open left half of the complex plane. Definition 2: The numerical scheme (1.4) is ( )A α − Stable for some 0, 2 πα  ∈    , if the wedge ( ){ }:| | , 0s z Arg z zα < α= − ≠ is contained in the region of absolute stability. The largest maxα is the angle of absolute stability. Definition 3: The numerical scheme (1.4) is stiffly stable if (i) it is absolutely stable in the Region { }1 :| Re(z) | LR z D= ≤ and (ii) accurate in the region ( ) ( ){ }2 1: D | Re | D ; | Im | D ,L RR z z z< < < = such American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 303 that the stability region is contained in the region 1 2R R∪ . The numerical scheme is Zero-Stable since the roots of the first characteristics polynomial ( ) 1k kr r rr −= − satisfy | r | 1i ≤ with roots of [ ] 1ir = being simple. To investigate the stability properties of the family of the hybrid multistep methods (1.4), we employ the boundary locus approach discussed in[14]. Substituting the hybrid predictors in (1.6) into (1.5) then into (1.4) at the hybrid points to yield a scheme, the resulting scheme for fixed k is applied to the scalar test problem ' ,y yλ= 2'' ,y yλ= ( )Re 0λ < which yields the stability polynomials as ( ) ( ) ( ) ( ) ( )( )1 2 0 , , m k m m mk k j k j p kv j r z r r z r H r z z rπ β β λ− =   = − − + −     ∑ (5.1) where ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )( )2 2 0 0 , ... ... ... l l k kl l l l l lk k j k p j j k kv v j j H r z r z r r z z r z Tβ β β β λ λ− − − = =       = − + + + +        ∑ ∑ and ( ) ( ) ( )2 0 k l l lj k j k k j T r z z rβ β λ− − − = = + +∑ The boundary plots are obtained from the stability polynomials for various k. 5. The Stability Plots of the hybrid method The following are the boundary plots of the implicit hybrid scheme derived in: The boundary loci reveal that the scheme (1.4) is zero-stable. For 6k ≤ , it is A -Stable and ( )A α -Stable for 6k > to k=9. American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 304 Figure 4 Figure 5 Figure 6 American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 305 6. Numerical implementations This section considers numerical implementation of the new hybrid methods (1.4) on some stiff initial value problems in ordinary differential equations. Since the method is an implicit method, the implicitness is resolved by applying the Newton scheme [ ] [ ] [ ]( ) [ ]( )11r r r r n k n k n k n ky y J y F y −+ + + + += − , 0,1,2,3,...r = (6.1) or a modification of (6.1) where [ ]( )r n kJ y + is the Jacobian matrix of the new hybrid method. The (6.1) requires starting value and is generated from the explicit scheme ( )1 1 , 2 2 r n n n n hy y f f p + −= + + = (6.2) Using fixed step-size h. The following problems are considered for implementation. Problem [1] The Chemical reaction problems in [17] 4 1 1 2 3' 0.04 10 ,y y y y = − + ( )1 0 1y = 4 7 2 2 1 2 3 2' 0.04 10 3.10 ,y y y y y= − − ( )2 0 0y = 7 2 3 2' 3.10 ,y y= ( )3 0 0y = 610h −= , [ ]0,3x∈ Problem [2] The non linear moderately stiff problems in [9] 1 1 2' 0.1 199.9 ,y y y= − − ( )1 0 2y = 2 2' 200y y= − , ( )2 0 1y = 0.0001h = with exact solution ( ) 0.1 200 1 x xy x e e− −= + and 200 2 (x) xy e−= American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 306 Problem [3] The Van der pol equation in [12] 1 2'y y= , ( )1 0 2y = ( )( )2 2 1 2 1' 1 /y y y y ε−= − , ( )2 0 0y = 0.001,h = 110ε −= 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 x-axis y- ax is ode15s1 ode15s2 ode15s3 y1 y2 y3 Figure 1: Graphical solution of problem1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.5 1 1.5 2 y1 y1E y2 y2E Figure 2: Graphical solution of problem 2 American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 307 0 5 10 15 -20 -15 -10 -5 0 5 10 15 20 x-axis fu nc tio n ax is data1 y1(ode15) y2 y2(ode15) Figure 3: Graphical solution of problem3 7. Conclusion This paper has presented a class of hybrid linear multistep methods (1.4) with nested hybrid predictors (1.5) for stiff initial value problems in ordinary differential equations. The hybrid scheme has high order stability and is seen to overcome Dahlquist order barrier on linear multistep methods (1.1). The scheme has been implemented on three stiff problems and the results in figures 1 and 3 show that the scheme (1.4) compares favourably with ODE15s of MATLAB in [13]. In figure 2, the graph is in alignment with the exact solution of the ODE. References [1] Brugnano, L & Trigiante, D.; Solving Differential Problems by Multistep Initial and Boundary Value Methods, Amsterdan: Gordon and Breach Science Publishers, 1998. [2] Butcher, J.C; A modified multistep method of numerical integration of ordinary Differential equations, J Ass. comput. Math; 1965, vol;12, pp.124-135. [3] Butcher, J.C; A Transformed implicit Runge-Kutta Method, J. Ass. comput. Math., 1979, Vol.26, pp.731-738. [4] Dahlquist, G. A; special stability problems for linear multistep methods, BIT. 1963, vol.3, pp. 27 [5] Donelson, J. & Hansen, E.; Cyclic Composite Multistep Predictor-Correctors Methods, SIAM, J. Num. Anal.,1971, Vol.8,pp.137-157. [6] Enright,W.H., Second Derivative Multistep Methods for stiff ODE’s, SAIM. J. Num.Anal., 1974, vol.11, ISS.2 pp. 321-331. [7] Enright, W.H., continuous numerical methods for ODE’s with defect control, J. computational. Appl. American Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2018) Volume 42, No 1, pp 297-308 308 math., Vol.25, (2000), pp. 159-170. [8] Esuabana I. M & Ekoro S. E.; Hybrid Linear Multistep Methods with Nested Hybrid Predictors for Solving Linear and Nonlinear Initial Value Problems in Ordinary Differential Equations, IISTE journal of Mathematical Theory and Modeling, 2017, vol. 7,iss. 11, pp. 77-88. [9] Fatunla, S. O.; Numerical Methods for Initial Value Problems in Ordinary Differential Equations, New York: Academic Press, 1988. [10] Gear, C. w.; Hybrid methods for IVP’s in ODEs, SIAM Journal on Numerical Analysis, vol.2,(1965), pp.69-86. [11] Gragg, W. B & Shetter, H. J.; Generalised Multistep Predictor-Correctors methods, J. Assoc. Comput. Mach.,1964, Vol.11, pp.188-209. [12] Hairer E. & Wanner G.; solving ordinary differential equation 11: Stiff and Differential Algebraic problems, 2nd rv.Ed. springer-verlag, New York,1996. [13] Higham, D. J.; Higham, N.J. –MATLAB Guide, Society of industrial and applied Mathematics (SIAM), Philadelphia, PA, 2000. [14] Ikhile, M. N. O & Okuonghae, R. I.; Stiffly Stable Continuous Extension of Second Derivative Linear Multistep Method with an off-step point for IVPs in ODEs, J. Nig. Assoc., Math. Phys., 2007, vol. 11, pp. 175-190. [15] Lambert, J. D.; Computational methods for Ordinary Differential Systems, Chichester; Wiley, 1973, pp.91. [16] Okuonghae, R. I., & Ikhile M. N. O., A class of Hybrid Linear Multistep Methods With ( )A α -Stable Properties for Stiff IVPs in ODEs, J. Num. Math. 2014, Vol. 8, iss. 4, pp. 441-469. [17] Robertson, H. H, The solution of a set of reaction rate equation in: Numerical Analysis: An introduction (J. Walsh, Ed.), academic Press, New York, 1966, pp. 178-182.