# Solving Nonlinear Fourth-Order Boundary Value Problems Using a Numerical Approach:(m+1)th-Step Block Method

Authors: Oluwaseun Adeyeye; Zurni Omar
Journal: International Journal of Differential Equations
 (2017)

---

## Abstract

Nonlinear boundary value problems (BVPs) are more tedious to solve than their linear counterparts. This is observed in the extra computation required when determining the missing conditions in transforming BVPs to initial value problems. Although a number of numerical approaches are already existent in literature to solve nonlinear BVPs, this article presents a new block method with improved accuracy to solve nonlinear BVPs. Am+1th-step block method is developed using a modified Taylor series approach to directly solve fourth-order nonlinear boundary value problems (BVPs) where m is the order of the differential equation under consideration. The schemes obtained were combined to simultaneously produce solution to the fourth-order nonlinear BVPs at m+1 points iteratively. The derived block method showed improved accuracy in comparison to previously existing authors when solving the same problems. In addition, the suitability of the m+1th-step block method was displayed in the solution for magnetohydrodynamic squeezing flow in porous medium.


---

## Body

## 1. Introduction

Boundary value problems (BVPs) arise in several branches of science ranging from physical sciences to engineering. There has been commendable progress in solving problems associated with nonlinear ordinary differential equations (ODEs) involving boundary conditions in recent years. These ODEs are sometimes needed to fulfil certain boundary conditions at more than one point of the independent variable which will result in the problem known as two-point boundary value problem. Two-point nonlinear BVPs often cannot be solved by analytical methods and thus finding approximate solutions for these problems becomes essential.This article considers the following special type of nonlinear boundary value problem:(1)yiv=fx,y,with the boundary conditions (2)ya=α1,y′a=α2,yb=β1,y′b=β2,where f is a continuous function on [a,b] and the parameters αi and βi for i=1,2 are constants.A variety of methods have been introduced to solve (1) such as shooting methods, splines methods, finite difference methods, finite element methods, differential transform methods, and collocation methods [1–4]. Recently, the adoption of various families of linear multistep method (LMM) for numerically approximating higher order ODEs has been proposed. However, some LMMs cannot directly solve these higher order ODEs and thus require reduction to a system of first-order ODEs. In some cases, the accuracy of LMMs is low such as the case of predictor-corrector methods which incur high computational rigour. This computational rigour involves the derivation of separate predictors for each grid point of the LMM as seen in the work of Kayode and Adeyeye [5] and Kayode and Obarhua [6]. These drawbacks caused for the introduction of block methods which were first proposed by Milne [7] as a means to obtain starting values for predictor-corrector methods. This concept was also further explored by Sarafyan [8].Block methods differ from alternate approaches such as differential transform method and collocation method. This is because the formulation of block methods is an evaluation of the linear multistep method at different grid points to generate a family of methods that can be applied to produce approximate solutions of ODEs at each grid point simultaneously. This advantage was mentioned by Lambert [9] among other advantages which include being self-starting, permitting easy change of step-length, and being less expensive in terms of function evaluations. Block methods also yield better accuracy when applied to numerical problems.This article introduces am+1th-step block method where m is the order of the differential equation. The block method is developed using a modification of the conventional Taylor series expansions approach by Lambert [9]. This derivation is shown in Section 2 of this article while Section 3 considers certain numerical examples and their results to show the accuracy of the block method.


## 2. Methodology

Lambert [9] highlighted three main approaches for developing LMMs. These include interpolation, numerical integration, and Taylor series expansions. This article adopts the Taylor series expansion approach for LMMs to develop the block method. However, certain modifications were introduced since Lambert [9] focused on first-order methods whereas this article develops a block method for fourth-order ODEs. Therefore, the approach was made suitable to develop block methods and not LMMs alone, hence the name modified Taylor series approach.


### 2.1. Derivation of them+1th-Step Block Method Using Modified Taylor Series Approach

The algorithm described below is used to show the steps involved in deriving them+1th-step block method using the modified Taylor series approach where m is the order of the differential equation.Algorithm 1.
Start.
Step 1. Obtain the coefficients of the initial multistep scheme: (3)yn+m+1=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,v=1,2,…,m.
Step 2. Obtain the coefficients of the additional schemes: (4)yn+jm+1=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,yn+jm+2=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,⋮yn+jk=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j.
Step 3. Derive the coefficients of the (1st,2nd,…,m-1th) derivative schemes: (5)ynϱ=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,yn+1ϱ=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,⋮yn+m+1ϱ=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,where ϱ=1,2,…,m-1.
Step 4. Combine schemes obtained in Steps 1, 2, and 3 above to form a system of equations with matrix form equivalent Ax=B where x=N0,N1,N2,…,Nm-1T and N0=yn+(m-3),yn+(m-2),…,yn+(m+1)T, N1=yn+(m-3)′,yn+(m-2)′,…,yn+(m+1)′T, N2=yn+(m-3)′′,yn+(m-2)′′,…,yn+(m+1)′′T, Nm-1=yn+(m-3)m-1,yn+(m-2)m-1,…,yn+(m+1)m-1T
Step 5. Adopt matrix inverse approach to system of equations in Step 4 to obtain the expected block method.
Stop.In Algorithm1, yn+a=y(xn+a)=y(xn+ah), fn+j=f(xn+j,yn+j), and αjv and βj are constants with v defined in Step 1 of Algorithm 1.Note that in Step1 of Algorithm 1, the expected αjv are αj1,αj2,…,αjm. The αjv-values can take kCv forms and jv-values not chosen will be used as evaluation points when developing the additional methods in Step 2.Steps1–3 of Algorithm 1 require expanding individual terms using Taylor Series expansion such as (6)yn=yxn,yn+1=yxn+h=yxn+hy′xn+h22!y′′xn+⋯,⋮yn+k-1=yxn+k-1h=yxn+k-1hy′xn+k-1h22!y′′xn+⋯,yn+k=yxn+kh=yxn+khy′xn+kh22!y′′xn+…,fn=ymxn,fn+1=ymxn+h=ymxn+hym+1xn+h22!ym+2xn+⋯,⋮fn+k-1=ymxn+k-1h=ymxn+k-1hym+1xn+k-1h22!ym+2xn+⋯,fn+k=ymxn+kh=ymxn+khym+1xn+kh22!ym+2xn+⋯.Substituting these expansions in individual equations and equating coefficients of y(m)xn presents the resulting expressions in matrix form (7)Ax=B,where (8)A=111…1000…00h2h…k-1h⋮⋮⋮…⋮0h22!2h22!…k-1h22!⋮⋮⋮…⋮⋮⋮⋮…⋮⋮⋮⋮…⋮0hmm!2hmm!…k-1hmm!111…10hm+1m+1!2hm+1m+1!…k-1hm+1m+1!0h2h…kh0..….0h22!2h22!…kh22!0..….0..….0..….0..….0..….0..….0h2k2k!2h2k2k!…k-1h2k2k!0h2k-m2k-m!2h2k-m2k-m!…kh2k-m2k-m!,x=α0,α1,α2,…,αk-1,β0,β1,β2,…,βkT,B=1,kh,kh22!,kh33!,kh44!,…,kh2k2k!T.Note that Algorithm 1 will not successfully obtain the required block method if matrix A is singular. Thus, the nonsingularity of the resulting matrices is discussed.


### 2.2. Nonsingularity of Resulting Matrices

The matrixA in (7) is a square matrix with det⁡(A)≠0 which follows from the theorems below.Theorem 2.
Suppose thatA is a square matrix with a row where every entry is zero, or a column where every entry is zero. Then det⁡(A)=0.Theorem 3.
Suppose thatA is a square matrix with two equal rows, or two equal columns. Then det⁡(A)=0.With respect to Theorem2, since the matrix A does not have a row or column where every entry is zero, then its inverse exists. On the other hand, matrix A has no equal rows or columns which further affirms that its inverse exists.Theorems2 and 3 are sufficient conditions to show that the inverse of the resulting matrix will always exist. In addition, the case of linear dependency is considered as defined in the following theorems.Theorem 4.
If matrixA has linearly dependent columns, then det⁡(A)=0.Theorem 5.
The rank of a matrixA equals the maximum number of linearly independent column vectors. The matrix A has the same number of linearly independent row vectors as it has linearly independent column vectorsThus, Theorem5 is tested for resulting matrices obtained in developing the m+1th-step block method to show that A-1 exists.


### 2.3. Specification of them+1th-Step Block Method

Following Algorithm1, the specification of the m+1th-step block method is as follows. From Step 1, the initial multistep scheme for the m+1th-step block method in terms of m is (9)yn+m+1=αm-4yn+m-4+αm-3yn+m-3+αm-2yn+m-2+αm-1yn+m-1+αmyn+m+βm-4fn+m-4+βm-3fn+m-3+βm-2fn+m-2+βm-1fn+m-1+βmfn+m+βm+1fn+m+1.Now, considering (9), the individual terms are expanded using Taylor series expansion as defined in (6). The resulting expansions are substituted back in (9) and rewritten in matrix form Ax=B, where (10)A=11111000000m-4hm-3hm-2hm-1hmh000000m-4h22!m-3h22!m-2h22!m-1h22!mh22!000000m-4h33!m-3h33!m-2h33!m-1h33!mh33!000000m-4h44!m-3h44!m-2h44!m-1h44!mh44!111111m-4h55!m-3h55!m-2h55!m-1h55!mh55!m-4hm-3hm-2hm-1hmhm+1hm-4h66!m-3h66!m-2h66!m-1h66!mh66!m-4h22!m-3h22!m-2h22!m-1h22!mh22!m+1h22!m-4h77!m-3h77!m-2h77!m-1h77!mh77!m-4h33!m-3h33!m-2h33!m-1h33!mh33!m+1h33!m-4h88!m-3h88!m-2h88!m-1h88!mh88!m-4h44!m-3h44!m-2h44!m-1h44!mh44!m+1h44!m-4h99!m-3h99!m-2h99!m-1h99!mh99!m-4h55!m-3h55!m-2h55!m-1h55!mh55!m+1h55!m-4h1010!m-3h1010!m-2h1010!m-1h1010!mh1010!m-4h66!m-3h66!m-2h66!m-1h66!mh66!m+1h66!,x=αm-4,αm-3,αm-2,αm-1,αm,βm-4,βm-3,βm-2,βm-1,βm,βm+1T,B=1,m+1h,m+1h22!,m+1h33!,m+1h44!,m+1h55!,m+1h66!,m+1h77!,m+1h88!,m+1h99!,m+1h1010!T,where matrix A has rank=11 which implies that there are no linearly dependent columns or rows and the inverse exists. This follows from the theorems in Section 2.2 showing that the matrix is nonsingular. Therefore, the scheme in (9) is obtained using matrix inverse method and substituting the value of m as (11)yn+m+1=yn+m-4-5yn+m-3+10yn+m-2-10yn+m-1+5yn+m+h4720fn+m-4-125fn+m-3-350fn+m-2+350fn+m-1+125fn+m+fn+m+1.Following the subsequent steps of Algorithm 1, the specification of the m+1th-step block method is as follows: (12)1814400yn+m-3=1814400yn+hyn′+h22yn′′+h36yn′′′+h449126fn+49045fn+m-3-40160fn+m-2+25430fn+m-1-9310fn+m+1469fn+m+1,14175yn+m-2=14175yn+2hyn′+2h2yn′′+4h33yn′′′+h44264fn+7960fn+m-3-4910fn+m-2+3080fn+m-1-1120fn+m+176fn+m+1,22400yn+m-1=22400yn+3hyn′+9h22yn′′+9h32yn′′′+h425488fn+63315fn+m-3-26460fn+m-2+19170fn+m-1-7020fn+m+1107fn+m+1,14175yn+m=14175yn+4hyn′+8h2yn′′+32h33yn′′′+h440448fn+116480fn+m-3-29440fn+m-2+33280fn+m-1-11360fn+m+1792fn+m+1,yn+m+1=yn+5hyn′+25h22yn′′+125h36yn′′′+h472576418250fn+1315625fn+m-3-175000fn+m-2+418750fn+m-1-106250fn+m+18625fn+m+1,(13)40320yn+m-3′=40320yn′+hyn′′+h22yn′′′+h33929fn+4975fn+m-3-3862fn+m-2+2422fn+m-1-883fn+m+139fn+m+1,630yn+m-2′=630yn′+2hyn′′+2h2yn′′′+h3317fn+734fn+m-3-380fn+m-2+244fn+m-1-89fn+m+14fn+m+1,4480yn+m-1′=4480yn′+3hyn′′+9h22yn′′′+h35481fn+16119fn+m-3-4374fn+m-2+4230fn+m-1-1539fn+m+243fn+m+1,315yn+m′=315yn′+4hyn′′+8h2yn′′′+h3712fn+2336fn+m-3-224fn+m-2+704fn+m-1-200fn+m+32fn+m+1,8064yn+m+1′=8064yn′+5hyn′′+25h22yn′′′+h329125fn+101875fn+m-3+1250fn+m-2+38750fn+m-1-4375fn+m+1375fn+m+1,(14)10080yn+m-3′′=10080yn′′+hyn′′′+h22462fn+4315fn+m-3-3044fn+m-2+1882fn+m-1-682fn+m+107fn+m+1,630yn+m-2′′=630yn′′+2hyn′′′+h2355fn+1088fn+m-3-370fn+m-2+272fn+m-1-101fn+m+16fn+m+1,1120yn+m-1′′=1120yn′′+3hyn′′′+h2984fn+3501fn+m-3-72fn+m-2+870fn+m-1-288fn+m+45fn+m+1,315yn+m′′=315yn′′+4hyn′′′+h2376fn+1424fn+m-3+176fn+m-2+608fn+m-1-80fn+m+16fn+m+1,2016yn+m+1′′=2016yn′′+5hyn′′′+h23050fn+11875fn+m-3+2500fn+m-2+6250fn+m-1+1250fn+m+275fn+m+1,(15)1440yn+1′′′=1440yn′′′+h475fn+1427fn+m-3-798fn+m-2+482fn+m-1-173fn+m+27fn+m+1,90yn+2′′′=90yn′′′+h28fn+129fn+m-3+14fn+m-2+14fn+m-1-6fn+m+fn+m+1,160yn+3′′′=160yn′′′+h51fn+219fn+m-3+114fn+m-2+114fn+m-1-21fn+m+3fn+m+1,45yn+4′′′=45yn′′′+h14fn+64fn+m-3+24fn+m-2+64fn+m-1+14fn+m,288yn+5′′′=288yn′′′+h95fn+375fn+m-3+250fn+m-2+250fn+m-1+375fn+m+95fn+m+1.


### 2.4. Order and Stability Properties of them+1th-Step Block Method

To ensure convergence of the block method, its consistency and zero-stability need to be investigated. This follows from Fatunla (1988) which states that a linear multistep method is convergent iff it is consistent and zero-stable.Starting with the consistency property, a linear multistep method isconsistent if it has order p≥1. Thus, the order of the m+1th-step block method is investigated.With reference to the definition in Lambert [9], Henrici [11], and Butcher [12], the order and error constant of the m+1th-step block method follow Definition 6.Definition 6.
The linear operator associated with LMM is defined as(16)Lyx;h=∑j=0kαjyn+j-βjfn+j.On expanding yn+j and fn+j to obtain (17)Lyx;h=C0yxn+C1hy′xn+⋯+Cqhqyqxn+⋯,where (18)C0=α0+α1+α2+⋯+αk,C1=α1+2α2+⋯+kαk,⋮Cq=1q!α1+2qα2+⋯+kqαk-1q-m!β1+2q-mβ2+⋯+kq-mβk,q=2,3,…,the method is said to be of order p if C0=C1=⋯=Cp=Cp+1=⋯=Cp+(m-1)=0, Cp+m≠0 and Cp+m is the error constant.The integrators of the block method (12) are of order six methods with the error constants, C10 obtained as 2323/3628800, 137/14175, 1737/44800, 1408/14175, and 29375/145152, respectively. Having order p>1, the consistency of the block method is affirmed.Moving on to the second criterion for convergence which is the zero-stability of the block method. Note that this is the most important stability property a good numerical method should possess as it ensures convergence. The key word “zero” is based on the stability phenomenon in terms of convergence in the limit as step-size (h) tends to zero.Therefore, to test the zero-stability of them+1th-step block method, the integrators are normalized to give the first characteristic polynomial ρ(r) as (19)ρr=det⁡rA0-A1=rmr-1,with A0=5×5 identity matrix(20)A1=0000100001000010000100001.The roots of ρ(r)=0 satisfy |rj|≤1. Hence, the m+1th-step block method is zero-stable.


## 2.1. Derivation of them+1th-Step Block Method Using Modified Taylor Series Approach

The algorithm described below is used to show the steps involved in deriving them+1th-step block method using the modified Taylor series approach where m is the order of the differential equation.Algorithm 1.
Start.
Step 1. Obtain the coefficients of the initial multistep scheme: (3)yn+m+1=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,v=1,2,…,m.
Step 2. Obtain the coefficients of the additional schemes: (4)yn+jm+1=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,yn+jm+2=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,⋮yn+jk=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j.
Step 3. Derive the coefficients of the (1st,2nd,…,m-1th) derivative schemes: (5)ynϱ=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,yn+1ϱ=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,⋮yn+m+1ϱ=∑j=m-4mαjvyn+j+∑j=m-4m+1βjfn+j,where ϱ=1,2,…,m-1.
Step 4. Combine schemes obtained in Steps 1, 2, and 3 above to form a system of equations with matrix form equivalent Ax=B where x=N0,N1,N2,…,Nm-1T and N0=yn+(m-3),yn+(m-2),…,yn+(m+1)T, N1=yn+(m-3)′,yn+(m-2)′,…,yn+(m+1)′T, N2=yn+(m-3)′′,yn+(m-2)′′,…,yn+(m+1)′′T, Nm-1=yn+(m-3)m-1,yn+(m-2)m-1,…,yn+(m+1)m-1T
Step 5. Adopt matrix inverse approach to system of equations in Step 4 to obtain the expected block method.
Stop.In Algorithm1, yn+a=y(xn+a)=y(xn+ah), fn+j=f(xn+j,yn+j), and αjv and βj are constants with v defined in Step 1 of Algorithm 1.Note that in Step1 of Algorithm 1, the expected αjv are αj1,αj2,…,αjm. The αjv-values can take kCv forms and jv-values not chosen will be used as evaluation points when developing the additional methods in Step 2.Steps1–3 of Algorithm 1 require expanding individual terms using Taylor Series expansion such as (6)yn=yxn,yn+1=yxn+h=yxn+hy′xn+h22!y′′xn+⋯,⋮yn+k-1=yxn+k-1h=yxn+k-1hy′xn+k-1h22!y′′xn+⋯,yn+k=yxn+kh=yxn+khy′xn+kh22!y′′xn+…,fn=ymxn,fn+1=ymxn+h=ymxn+hym+1xn+h22!ym+2xn+⋯,⋮fn+k-1=ymxn+k-1h=ymxn+k-1hym+1xn+k-1h22!ym+2xn+⋯,fn+k=ymxn+kh=ymxn+khym+1xn+kh22!ym+2xn+⋯.Substituting these expansions in individual equations and equating coefficients of y(m)xn presents the resulting expressions in matrix form (7)Ax=B,where (8)A=111…1000…00h2h…k-1h⋮⋮⋮…⋮0h22!2h22!…k-1h22!⋮⋮⋮…⋮⋮⋮⋮…⋮⋮⋮⋮…⋮0hmm!2hmm!…k-1hmm!111…10hm+1m+1!2hm+1m+1!…k-1hm+1m+1!0h2h…kh0..….0h22!2h22!…kh22!0..….0..….0..….0..….0..….0..….0h2k2k!2h2k2k!…k-1h2k2k!0h2k-m2k-m!2h2k-m2k-m!…kh2k-m2k-m!,x=α0,α1,α2,…,αk-1,β0,β1,β2,…,βkT,B=1,kh,kh22!,kh33!,kh44!,…,kh2k2k!T.Note that Algorithm 1 will not successfully obtain the required block method if matrix A is singular. Thus, the nonsingularity of the resulting matrices is discussed.


## 2.2. Nonsingularity of Resulting Matrices

The matrixA in (7) is a square matrix with det⁡(A)≠0 which follows from the theorems below.Theorem 2.
Suppose thatA is a square matrix with a row where every entry is zero, or a column where every entry is zero. Then det⁡(A)=0.Theorem 3.
Suppose thatA is a square matrix with two equal rows, or two equal columns. Then det⁡(A)=0.With respect to Theorem2, since the matrix A does not have a row or column where every entry is zero, then its inverse exists. On the other hand, matrix A has no equal rows or columns which further affirms that its inverse exists.Theorems2 and 3 are sufficient conditions to show that the inverse of the resulting matrix will always exist. In addition, the case of linear dependency is considered as defined in the following theorems.Theorem 4.
If matrixA has linearly dependent columns, then det⁡(A)=0.Theorem 5.
The rank of a matrixA equals the maximum number of linearly independent column vectors. The matrix A has the same number of linearly independent row vectors as it has linearly independent column vectorsThus, Theorem5 is tested for resulting matrices obtained in developing the m+1th-step block method to show that A-1 exists.


## 2.3. Specification of them+1th-Step Block Method

Following Algorithm1, the specification of the m+1th-step block method is as follows. From Step 1, the initial multistep scheme for the m+1th-step block method in terms of m is (9)yn+m+1=αm-4yn+m-4+αm-3yn+m-3+αm-2yn+m-2+αm-1yn+m-1+αmyn+m+βm-4fn+m-4+βm-3fn+m-3+βm-2fn+m-2+βm-1fn+m-1+βmfn+m+βm+1fn+m+1.Now, considering (9), the individual terms are expanded using Taylor series expansion as defined in (6). The resulting expansions are substituted back in (9) and rewritten in matrix form Ax=B, where (10)A=11111000000m-4hm-3hm-2hm-1hmh000000m-4h22!m-3h22!m-2h22!m-1h22!mh22!000000m-4h33!m-3h33!m-2h33!m-1h33!mh33!000000m-4h44!m-3h44!m-2h44!m-1h44!mh44!111111m-4h55!m-3h55!m-2h55!m-1h55!mh55!m-4hm-3hm-2hm-1hmhm+1hm-4h66!m-3h66!m-2h66!m-1h66!mh66!m-4h22!m-3h22!m-2h22!m-1h22!mh22!m+1h22!m-4h77!m-3h77!m-2h77!m-1h77!mh77!m-4h33!m-3h33!m-2h33!m-1h33!mh33!m+1h33!m-4h88!m-3h88!m-2h88!m-1h88!mh88!m-4h44!m-3h44!m-2h44!m-1h44!mh44!m+1h44!m-4h99!m-3h99!m-2h99!m-1h99!mh99!m-4h55!m-3h55!m-2h55!m-1h55!mh55!m+1h55!m-4h1010!m-3h1010!m-2h1010!m-1h1010!mh1010!m-4h66!m-3h66!m-2h66!m-1h66!mh66!m+1h66!,x=αm-4,αm-3,αm-2,αm-1,αm,βm-4,βm-3,βm-2,βm-1,βm,βm+1T,B=1,m+1h,m+1h22!,m+1h33!,m+1h44!,m+1h55!,m+1h66!,m+1h77!,m+1h88!,m+1h99!,m+1h1010!T,where matrix A has rank=11 which implies that there are no linearly dependent columns or rows and the inverse exists. This follows from the theorems in Section 2.2 showing that the matrix is nonsingular. Therefore, the scheme in (9) is obtained using matrix inverse method and substituting the value of m as (11)yn+m+1=yn+m-4-5yn+m-3+10yn+m-2-10yn+m-1+5yn+m+h4720fn+m-4-125fn+m-3-350fn+m-2+350fn+m-1+125fn+m+fn+m+1.Following the subsequent steps of Algorithm 1, the specification of the m+1th-step block method is as follows: (12)1814400yn+m-3=1814400yn+hyn′+h22yn′′+h36yn′′′+h449126fn+49045fn+m-3-40160fn+m-2+25430fn+m-1-9310fn+m+1469fn+m+1,14175yn+m-2=14175yn+2hyn′+2h2yn′′+4h33yn′′′+h44264fn+7960fn+m-3-4910fn+m-2+3080fn+m-1-1120fn+m+176fn+m+1,22400yn+m-1=22400yn+3hyn′+9h22yn′′+9h32yn′′′+h425488fn+63315fn+m-3-26460fn+m-2+19170fn+m-1-7020fn+m+1107fn+m+1,14175yn+m=14175yn+4hyn′+8h2yn′′+32h33yn′′′+h440448fn+116480fn+m-3-29440fn+m-2+33280fn+m-1-11360fn+m+1792fn+m+1,yn+m+1=yn+5hyn′+25h22yn′′+125h36yn′′′+h472576418250fn+1315625fn+m-3-175000fn+m-2+418750fn+m-1-106250fn+m+18625fn+m+1,(13)40320yn+m-3′=40320yn′+hyn′′+h22yn′′′+h33929fn+4975fn+m-3-3862fn+m-2+2422fn+m-1-883fn+m+139fn+m+1,630yn+m-2′=630yn′+2hyn′′+2h2yn′′′+h3317fn+734fn+m-3-380fn+m-2+244fn+m-1-89fn+m+14fn+m+1,4480yn+m-1′=4480yn′+3hyn′′+9h22yn′′′+h35481fn+16119fn+m-3-4374fn+m-2+4230fn+m-1-1539fn+m+243fn+m+1,315yn+m′=315yn′+4hyn′′+8h2yn′′′+h3712fn+2336fn+m-3-224fn+m-2+704fn+m-1-200fn+m+32fn+m+1,8064yn+m+1′=8064yn′+5hyn′′+25h22yn′′′+h329125fn+101875fn+m-3+1250fn+m-2+38750fn+m-1-4375fn+m+1375fn+m+1,(14)10080yn+m-3′′=10080yn′′+hyn′′′+h22462fn+4315fn+m-3-3044fn+m-2+1882fn+m-1-682fn+m+107fn+m+1,630yn+m-2′′=630yn′′+2hyn′′′+h2355fn+1088fn+m-3-370fn+m-2+272fn+m-1-101fn+m+16fn+m+1,1120yn+m-1′′=1120yn′′+3hyn′′′+h2984fn+3501fn+m-3-72fn+m-2+870fn+m-1-288fn+m+45fn+m+1,315yn+m′′=315yn′′+4hyn′′′+h2376fn+1424fn+m-3+176fn+m-2+608fn+m-1-80fn+m+16fn+m+1,2016yn+m+1′′=2016yn′′+5hyn′′′+h23050fn+11875fn+m-3+2500fn+m-2+6250fn+m-1+1250fn+m+275fn+m+1,(15)1440yn+1′′′=1440yn′′′+h475fn+1427fn+m-3-798fn+m-2+482fn+m-1-173fn+m+27fn+m+1,90yn+2′′′=90yn′′′+h28fn+129fn+m-3+14fn+m-2+14fn+m-1-6fn+m+fn+m+1,160yn+3′′′=160yn′′′+h51fn+219fn+m-3+114fn+m-2+114fn+m-1-21fn+m+3fn+m+1,45yn+4′′′=45yn′′′+h14fn+64fn+m-3+24fn+m-2+64fn+m-1+14fn+m,288yn+5′′′=288yn′′′+h95fn+375fn+m-3+250fn+m-2+250fn+m-1+375fn+m+95fn+m+1.


## 2.4. Order and Stability Properties of them+1th-Step Block Method

To ensure convergence of the block method, its consistency and zero-stability need to be investigated. This follows from Fatunla (1988) which states that a linear multistep method is convergent iff it is consistent and zero-stable.Starting with the consistency property, a linear multistep method isconsistent if it has order p≥1. Thus, the order of the m+1th-step block method is investigated.With reference to the definition in Lambert [9], Henrici [11], and Butcher [12], the order and error constant of the m+1th-step block method follow Definition 6.Definition 6.
The linear operator associated with LMM is defined as(16)Lyx;h=∑j=0kαjyn+j-βjfn+j.On expanding yn+j and fn+j to obtain (17)Lyx;h=C0yxn+C1hy′xn+⋯+Cqhqyqxn+⋯,where (18)C0=α0+α1+α2+⋯+αk,C1=α1+2α2+⋯+kαk,⋮Cq=1q!α1+2qα2+⋯+kqαk-1q-m!β1+2q-mβ2+⋯+kq-mβk,q=2,3,…,the method is said to be of order p if C0=C1=⋯=Cp=Cp+1=⋯=Cp+(m-1)=0, Cp+m≠0 and Cp+m is the error constant.The integrators of the block method (12) are of order six methods with the error constants, C10 obtained as 2323/3628800, 137/14175, 1737/44800, 1408/14175, and 29375/145152, respectively. Having order p>1, the consistency of the block method is affirmed.Moving on to the second criterion for convergence which is the zero-stability of the block method. Note that this is the most important stability property a good numerical method should possess as it ensures convergence. The key word “zero” is based on the stability phenomenon in terms of convergence in the limit as step-size (h) tends to zero.Therefore, to test the zero-stability of them+1th-step block method, the integrators are normalized to give the first characteristic polynomial ρ(r) as (19)ρr=det⁡rA0-A1=rmr-1,with A0=5×5 identity matrix(20)A1=0000100001000010000100001.The roots of ρ(r)=0 satisfy |rj|≤1. Hence, the m+1th-step block method is zero-stable.


## 3. Results and Discussion

This section tests them+1th-step block method on some nonlinear problems. The numerical results are shown in Tables 1–3 and Figures 1–3.Table 1
Comparison of the(m+1)th-step block method with Mustafa et al. [10] for solving Example 7.
x
Exact solution
Computed solution
Error [10]
Error ((m+1)th-step block method)
0.0
0.00000000000
0.00000000000
0.000000e+00
0.000000e+00
0.1
0.09531017980
0.09531018728
0.0002954265
7.472134e-09
0.2
0.18232155679
0.18232159067
0.0008719341
3.387990e-08
0.3
0.26236426447
0.26236434273
0.0014096072
7.826736e-08
0.4
0.33647223662
0.33647237832
0.0017352146
1.417001e-07
0.5
0.40546510811
0.40546533420
0.0017810699
2.260955e-07
0.6
0.47000362925
0.47000394905
0.0015577013
3.198013e-07
0.7
0.53062825106
0.53062863602
0.0011349902
3.849590e-07
0.8
0.58778666490
0.58778704558
0.0006286279
3.806745e-07
0.9
0.64185388617
0.64185415218
0.0001902154
2.660067e-07
1.0
0.69314718056
0.69314718056
0.000000e+00
0.000000e+00Table 2
Comparison of the exact and computed solution of Example8.
x
Exact solution
Computed solution
Error [10]
Error ((m+1)th-step block method)
0.0
0.00000000000
0.00000000000
0.000000e+00
0.000000e+00
0.1
0.01981000000
0.01981000000
0.0004095
0.000000e+00
0.2
0.07712000000
0.07712000000
0.0025752
0.000000e+00
0.3
0.16623000000
0.16623000000
0.0066432
0.000000e+00
0.4
0.27904000000
0.27904000000
0.0115595
0.000000e+00
0.5
0.40625000000
0.40625000000
0.0156708
0.000000e+00
0.6
0.53856000000
0.53856000000
0.0173246
0.000000e+00
0.7
0.66787000000
0.66787000000
0.0154706
0.000000e+00
0.8
0.78848000000
0.78848000000
0.0102612
0.000000e+00
0.9
0.89829000000
0.89829000000
0.0036517
0.000000e+00
1.0
1.00000000000
1.00000000000
0.000000e+00
0.000000e+00Table 3
Absolute errors of fifth-order HAM and(m+1)th-step block method when Rmp=0.
z
Exact solution
Error [13]
Error ((m+1)th-step block method)
0.0
0
0
0
0.1
0.085233703438701791
7.58785506649317×10-10
4.032885×10-14
0.2
0.171320454429454980
1.39478356642186×10-9
2.363387×10-13
0.3
0.259121838110931650
1.80948145356296×10-9
6.848411×10-13
0.4
0.349516600242079760
1.94815263920844×10-9
1.489975×10-12
0.5
0.443409441985037010
1.81075676675135×10-9
2.768896×10-12
0.6
0.541740074458440520
1.45204859247627×10-9
4.502620×10-12
0.7
0.645492623682151550
9.70818092582703×10-10
6.029954×10-12
0.8
0.755705480041236500
4.90335771985428×10-10
6.408873×10-12
0.9
0.873481690845957730
1.33847599670389×10-10
4.708123×10-12
1.0
1.000000000000000000
2.22044604925031×10-16
0.000000×100Figure 1
Comparison of the exact and computed solution of Example7.Figure 2
Comparison of the exact and computed solution of Example8.Figure 3
Comparison of solution between HAM [13] and m+1th-step block method.Example 7.
Consider the following nonlinear boundary value problem [10]: (21)yiv-6exp⁡-4y=-121+x4,with boundary conditions (22)y0=0,y′0=1,y1=ln⁡2=y′1=0.5.The exact solution of Example 7 is y=ln⁡(1+x). The obtained numerical results for this problem are presented in Table 1 with h=10-1. The maximum absolute error obtained by the m+1th-step block method is 3.84959×10-7 which is more accurate than the maximum error of 1.78×10-3 by Mustafa et al. [10]. The graphical comparison between exact and computed solution is shown in Figure 1.Example 8.
Consider the following nonlinear boundary value problem [10]: (23)yiv=y2-x10+-4x8-4x7+8x6-4x4+120x-48,with boundary conditions (24)y0=y′0=0,y1=y′1=1.The exact solution of Example 8 is y=x5-2x4+2x2. The obtained numerical results for this problem are presented in Table 2 with h=10-1. The m+1th-step block method gives precise and accurate results as the exact solution. This is far more encouraging than the maximum error of 1.73×10-2 by Mustafa et al. [10]. The graphical comparison between exact and computed solution is also shown in Figure 2.Examples7 and 8 considered the solution of nonlinear boundary value problems solved by Mustafa et al. [10]. In their work, the authors adopted a numerical approach based on subdivision schemes. Although their approach gave good results, the m+1th-step block method gave better results in terms of accuracy. This superiority in accuracy of the m+1th-step block method is resultant from its self-starting implementation approach instead of the approach of Mustafa et al. [10] requiring choosing different subdivision schemes with certain adjustment of boundary conditions. The self-starting approach of the block method requires no starting values which could reduce the accuracy of the method while also increasing the computational rigour. Rather, the integrators of m+1th-step block method were combined as direct simultaneous integrators for the solution of the nonlinear boundary value problems.In addition, to further show the suitability of them+1th-step block method, a physical problem is solved and the results are compared to exiting solutions in literature.


### 3.1.m+1th-Step Block Method Solution for Magnetohydrodynamic Squeezing Flow in Porous Medium

The study of squeezing effect, in addition to other properties such as magnetohydrodynamics (MHD) and porosity, has become one of the most active topics in fluid mechanics. Ullah et al. [13] made an effort to investigate MHD squeezing flow of Newtonian fluid between two parallel plates passing through porous medium by homotopy analysis method (HAM). The authors used similarity transforms to convert the governing partial differential equations to equivalent nonlinear boundary value problems and then solved by HAM. This article takes a step further to solve the resulting nonlinear boundary value problems using the m+1th-step block method.The resulting differential equation is a fourth-order nonlinear boundary value problem of the form(25)d4dz4fz+Rmpfzd3dz3fz-mpd2dz2fz-mhd2dz2fz=0,with boundary conditions (26)d2dz2f0=0,ddzf1=0,f0=0,f1=1,where Rmp is Reynold number and mh, mp are Hartmann numbers.Due to the difficulty to compute an exact solution for (25), Ullah et al. [13] computed varying solutions of (25) for different HAM orders. The m+1th-step block method is likewise adopted to solve (25) and convergence is observed to the results proposed by Ullah et al. [13] as seen in Figure 3.In addition, a special case of (25) is studied when the Reynold number is zero with exact solution obtained using the boundary conditions (26) as (27)fz=e2-MxeM-eM+Mxe-Mx+Mxe-Mxe2Me-MxMe2M-e2M+M+1,where M=mh+mp.Comparison is made between the exact solution, fifth-order HAM solution [13], and the m+1th-step block method in terms of absolute error.From Table3, an improved accuracy was displayed by the m+1th-step block method over the homotopy analysis method. This shows the block method is appropriate to evaluate the numerical solution of physical problems modelled as fourth-order nonlinear boundary value problems.


## 3.1.m+1th-Step Block Method Solution for Magnetohydrodynamic Squeezing Flow in Porous Medium

The study of squeezing effect, in addition to other properties such as magnetohydrodynamics (MHD) and porosity, has become one of the most active topics in fluid mechanics. Ullah et al. [13] made an effort to investigate MHD squeezing flow of Newtonian fluid between two parallel plates passing through porous medium by homotopy analysis method (HAM). The authors used similarity transforms to convert the governing partial differential equations to equivalent nonlinear boundary value problems and then solved by HAM. This article takes a step further to solve the resulting nonlinear boundary value problems using the m+1th-step block method.The resulting differential equation is a fourth-order nonlinear boundary value problem of the form(25)d4dz4fz+Rmpfzd3dz3fz-mpd2dz2fz-mhd2dz2fz=0,with boundary conditions (26)d2dz2f0=0,ddzf1=0,f0=0,f1=1,where Rmp is Reynold number and mh, mp are Hartmann numbers.Due to the difficulty to compute an exact solution for (25), Ullah et al. [13] computed varying solutions of (25) for different HAM orders. The m+1th-step block method is likewise adopted to solve (25) and convergence is observed to the results proposed by Ullah et al. [13] as seen in Figure 3.In addition, a special case of (25) is studied when the Reynold number is zero with exact solution obtained using the boundary conditions (26) as (27)fz=e2-MxeM-eM+Mxe-Mx+Mxe-Mxe2Me-MxMe2M-e2M+M+1,where M=mh+mp.Comparison is made between the exact solution, fifth-order HAM solution [13], and the m+1th-step block method in terms of absolute error.From Table3, an improved accuracy was displayed by the m+1th-step block method over the homotopy analysis method. This shows the block method is appropriate to evaluate the numerical solution of physical problems modelled as fourth-order nonlinear boundary value problems.


## 4. Conclusion

This article has introduced a numerical approach based on block methods derived using modified Taylor series approach. Them+1th-step block method was adopted for the numerical solution of different nonlinear fourth-order boundary value problems. The numerical results show that the impressive accuracy of the m+1th-step block method having the same computations as the exact solution is obtained as shown in Tables 1–3 and Figures 1–3. This grounds the suitability of the m+1th-step block method for solving fourth-order nonlinear boundary value problems. In addition, the suitability of the block method in application to physical problems was also investigated by presenting a solution to MHD squeezing flow in a porous medium. Convergence in solution and improved accuracy were properties displayed by the m+1th-step block method when solving this fluid model. Thus, the m+1th-step block method is appropriate for solving fourth-order nonlinear boundary value problems.



---

*Source: 4925914-2017-10-29.xml*