Higher Order Compact Finite Difference Schemes for Unsteady Boundary Layer Flow Problems

Abstract

We investigate the applicability of the compact finite difference relaxation method (CFDRM) in solving unsteady boundary layer flow problems modelled by nonlinear partial differential equations. The CFDRM utilizes the Gauss-Seidel approach of decoupling algebraic equations to linearize the governing equations and solve the resulting system of ordinary differential equations using compact finite difference schemes. The CFDRM has only been used to solve ordinary differential equations modelling boundary layer problems. This work extends its applications to nonlinear partial differential equations modelling unsteady boundary layer flows. The CFDRM is validated on two examples and the results are compared to results of the Keller-box method.

(ProQuest: ... denotes non-US-ASCII text omitted.)

P. G. Dlamini 1 and S. S. Motsa 2 and M. Khumalo 1

Academic Editor:Waqar Khan

1, Department of Mathematics, University of Johannesburg, P.O. Box 17011, Doornfontein 2028, South Africa
2, School of Mathematical Sciences, University of KwaZulu-Natal, Private Bag X01, Pietermaritzburg, Scottsville 3209, South Africa

Received 1 October 2013; Accepted 2 December 2013

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

The many, important applications associated with boundary layer flow and heat transfer induced by a stretching surface have made them one of the most studied problems in the field of fluid dynamics. These applications include aerodynamic extrusion of plastic sheets, the cooling process of metallic plates in a cooling bath, and the cooling and/or drying of paper and textiles Sakiadis [1] was the first to study boundary layer flows induced by a stretching surface. Since then, due to their industrial importance they have been studied by quite a number of researchers. Initially, steady flows were considered, but since it is not always to maintain steady state conditions, in recent years the focus has shifted towards unsteady flows. In this work, we consider unsteady flows modelled by partial differential equations.

In solving the unsteady boundary layer problems many researchers have employed the Williams III and Rhyne [2] transformation. Williams III and Rhyne introduced self-similar coordinates to transform the equations governing the unsteady flows by converting the infinite time scale of unsteady problems to a finite region of integration. Various analytical and numerical methods have been used to solve the transformed equations. The Keller-box method of Cebeci and Bradshaw [3] has been very popular in solving this kind of problems amongst researchers. These include Seshadri et al. [4] who combined the Keller-box and perturbation series approach for the solution of unsteady mixed convection flow along a heated vertical plate. Nazar et al. [5, 6] solved the unsteady boundary layer flow problem due to an impulsively stretching surface in a rotating fluid by means of the Keller-box numerical method, and they obtained a first-order perturbation approximation of the solution. The limitation of the perturbation approach, as noted by Liao [7], is that it gives solutions that are only valid for small time. Liao proposed the use of the homotopy analysis method (HAM) instead and since then the HAM has been dominantly used in solving the unsteady boundary layer flow problems (see [8-16]). The advantage of the HAM over the perturbation techniques is that it is able to produce solutions valid for all time.

In this work we propose the use of a new method that combines the use of higher order compact finite difference (CFD) schemes and the Gauss-Seidel approach of solving algebraic equations. The method called compact finite difference relaxation method (CFDRM) has only been used once by Dlamini et al. [17] to solve ordinary differential equations. The CFDRM is based on simple decoupling and rearrangement of the governing equations and numerically integrating the resulting equations using the CFD schemes. The use of higher order CFD schemes is motivated by the fact that they produce highly accurate solutions on coarser grids with greater computational efficiency. Various CFD schemes used for applications such as interpolation, filtering, and evaluating high-order derivatives were discussed in detail by Lele [18]. CFD schemes have been used to solve, for example, Burger's equation [19, 20], Navier-Stokes equation [21], Korteweg-de Vries equation [22], Black-Scholes equation [23], and many more [24-26]. In this work we explore their usage in solving unsteady boundary layer flow problems.

We consider two examples, the unsteady one-dimensional MHD laminar boundary layer flow due to an impulsively stretching surface that was previously studied by Srinivasa and Eswara [27], and the unsteady three-dimensional MHD flow and heat transfer over an impulsively stretching plate previously studied by Xu et al. [8]. The main objective of the study is to investigate the applicability of the CFDRM to solve unsteady boundary layer problems modelled by partial differential equations. The results are compared with results of the Keller-box method.

The rest of the paper is organized as follows. In Section 2, we discuss the derivation of the compact finite difference schemes, In Section 3 we discuss the development of the CFDRM for computing the solution of an unsteady MHD laminar boundary layer flow due to an impulsively stretching surface. Section 4 presents CFDRM implementation of an unsteady three-dimensional MHD flow and heat transfer over an impulsively stretching plate. Section 5 contains the results and discussion and the conclusions are given in Section 6.

2. Compact Finite Difference Schemes

In the derivation of the CFD schemes we consider a one-dimensional uniform mesh on the region [a,b] with nodes _xi (i=1,2,...,N) , where [figure omitted; refer to PDF] and a corresponding function _yi =y(_xi ) at the nodes. The distance between any two successive nodes is a constant h=_xi -_xi-1 . Sixth-order approximations of the first, second, and third derivatives at interior nodes can be obtained using the following schemes (see [18] for details): [figure omitted; refer to PDF]

For illustrative purposes we describe the application of the CFD schemes to second-order differential equations for y(x) with known boundary conditions at y(a) and y(b) . Consider the nonlinear differential equations [figure omitted; refer to PDF] where f(x,y,^{y[variant prime]} ) is a nonlinear function, p(x) and q(x) are known functions of x , and _ya and _yb are known constants. In solving (4) we apply the CFD approximation for the first and second derivatives given by (2) and (3), respectively, at the interior nodes (i=2,...,N-1) . Since we know boundary conditions at i=1 and i=N , the CFD schemes must be adjusted for the nodes near the boundary points. In order to maintain the order O(^h6 ) accuracy at the boundary points as in the interior points and to maintain the same tridiagonal format, we use the following one sided scheme at i=2 : [figure omitted; refer to PDF] and when i=N-1 we use [figure omitted; refer to PDF] where _ai ,_bi (i=1,...,7) are constants to be determined. To obtain a sixth-order accurate scheme we use Taylor series expansion about _x1 and _xN up to O(^h7 ) , in (5) and (6), respectively, and equate terms of order h . In each case, we obtain a system of seven linear algebraic equations in seven unknowns which are solved to give [figure omitted; refer to PDF]

Similarly, for the second derivatives, we use [figure omitted; refer to PDF] at i=2 and [figure omitted; refer to PDF] at i=N-1 . The parameters _ci , _di can be determined by expanding (8) and (9) using Taylor series and equating powers of h and subsequently solving the resulting equations. This gives [figure omitted; refer to PDF]

Using the above equations, the equations for approximating the first- and second-order derivatives can be expressed as [figure omitted; refer to PDF] where [figure omitted; refer to PDF]

From (11), ^{Y[variant prime]} and ^{Y[variant prime][variant prime]} can be expressed as [figure omitted; refer to PDF] where [figure omitted; refer to PDF]

3. Unsteady MHD Laminar Boundary Layer Flow Due to an Impulsively Stretching Surface

We consider the unsteady, laminar incompressible flow of a vicious electrically conducting fluid over a linear surface. This problem was studied by Srinivasa and Eswara [27]. They used the Keller-box method to solve the problem. After applying Williams III and Rhyne [2] similarity transformation to the original equation governing this flow they obtained the following dimensionless governing equation for the flow [figure omitted; refer to PDF] subject to the boundary conditions [figure omitted; refer to PDF] where ξ∈[0,1] is the dimensionless time scale defined as [figure omitted; refer to PDF] with b as positive constant and t as the time variable. In the analysis of boundary layer flow problems, a quantity that is of physical interest is the skin friction which in this model is given [4, 5, 7], in dimensionless form, as [figure omitted; refer to PDF] where R_ex is the local Reynolds number.

The initial unsteady solution at ξ=0 (τ=0 ) for the governing equation (15) is obtained as a solution of the equation [figure omitted; refer to PDF] where the primes denote differentiation with respect to η . Solving (19) gives [figure omitted; refer to PDF] where erfc(η) is the standard complementary error function defined by [figure omitted; refer to PDF]

The steady state solution when ξ=1 , corresponding to τ[arrow right]+∞ , is obtained from [figure omitted; refer to PDF] The solution to the above equation is [figure omitted; refer to PDF]

Now we discuss the development of the CFDRM to solve the partial differential equation (15). The CFDRM works for systems of differential equations and hence it is convenient to reduce the order of (15) from three to two. We set ^{f[variant prime]} =u , so that (15) becomes [figure omitted; refer to PDF]

The CFDRM algorithm uses the same ideas as the SRM described in [28-30]. Both employ the Gauss Seidel idea of decoupling systems of algebraic equations to decouple differential equations. In the framework of the SRM the resulting linear system of differential equations is integrated by the spectral collocation method. Instead of using spectral method, the CFDRM uses the compact finite differences to integrate the linearized system of differential equations. Applying the CFDRM on (24) we obtain the following linear PDEs: [figure omitted; refer to PDF] where [figure omitted; refer to PDF]

Here, r denotes the iteration level. All linear terms are evaluated at current iteration level (r+1 ) and nonlinear terms at the previous iteration level (r ). The initial approximations for solving (25)-(26) are obtained as the solutions at ξ=0 . Thus _f0 (η,ξ) and _u0 (η,ξ) are given by [figure omitted; refer to PDF] Starting from given initial approximations (29), the iteration scheme (25) can be solved iteratively for _ur+1 (η,ξ) when r=0,1,2,... . The solution for _ur+1 is used in (26) which is, in turn, solved for _fr+1 . To solve (25) we use sixth order compact finite difference schemes described in the previous section and use an implicit finite difference method in the ξ -direction.

The finite difference scheme with centering about a midpoint halfway between ^ξn+1 and ^ξn is applied in the η direction. This midpoint is defined as ^ξn+(1/2) =(^ξn+1 +^ξn )/2 . Applying the centering about ^ξn+(1/2) to any function, say u(η,ξ) and its associated derivative, we obtain [figure omitted; refer to PDF]

We first apply the CFD schemes on (25) and (26) to obtain [figure omitted; refer to PDF] with initial approximations [figure omitted; refer to PDF]

Next, we apply the finite difference scheme on (31) in the ξ -direction with centering about the midpoint ^ξn+(1/2) to obtain [figure omitted; refer to PDF] subject to the following boundary and initial conditions: [figure omitted; refer to PDF] where [figure omitted; refer to PDF] where I is an (_Nx +1)×(_Nx +1) identity matrix. The matrices _E1 ,_E2 ,_H1 , and _H2 are as described in Section 3 above. Thus, starting from the initial condition ^Ur+10 , given by (29), (35) can be solved iteratively to give approximate solutions for _ur+1 (η,ξ) , r=0,1,2,... , until a solution that converges to within a given accuracy level is obtained. The solution _ur+1 is used in (32) which is, in turn, solved for _fr+1 .

4. Unsteady Three-Dimensional MHD Flow and Heat Transfer over an Impulsively Stretching Plate

We consider the unsteady, three-dimensional flow of an electrically conducting fluid caused by an impulsive stretching surface in two lateral directions. The fluid is electrically conducting in a traverse magnetic field. The magnetic field is neglected under the assumption of a small magnetic Reynolds number. Xu et al. [8] solved this problem using the homotopy analysis method. Here we use the compact finite difference relaxation method to solve the governing equations. The flow is governed by the following three partial differential equations: [figure omitted; refer to PDF] with the following boundary conditions: [figure omitted; refer to PDF]

In the above equations prime denotes the derivative with respect to η and c the stretching parameter is a positive constant. M is the local Hartman number and Pr the Prandtl number. The initial unsteady solution can be found exactly by setting ξ=0 in the above equations and solving the resulting equations. The closed form analytical solutions are given by [figure omitted; refer to PDF]

Now we discuss the development of the spectral relaxation method to solve the system of partial differential equations (38)-(40). First, we set ^{f[variant prime]} =u and ^{s[variant prime]} =v , so that (38) and (39) become [figure omitted; refer to PDF] Applying the CFDRM on the resulting system of nonlinear partial differential equations gives the following linear partial differential equations: [figure omitted; refer to PDF] with the following boundary conditions: [figure omitted; refer to PDF] where [figure omitted; refer to PDF]

Starting from given initial approximations, denoted by _u0 (η,ξ) , _f0 (η,ξ) , _v0 (η,ξ) , _s0 (η,ξ) , and _g0 (η,ξ) , (44) can be solved iteratively for the unknown functions. To solve the above decoupled system of differential equations we apply CFD schemes on the space variable and finite differences in the time variable as described previously and obtain the following system of decoupled equations: [figure omitted; refer to PDF] with the following boundary conditions: [figure omitted; refer to PDF] and initial functions [figure omitted; refer to PDF] The matrices above are defined as [figure omitted; refer to PDF]

The vectors U , F , V , S , and G are the vectors of the functions u , f , v , s , and g , respectively, when evaluated at the grid points _ηj (j=0,1,...,_Nx ) .

5. Results and Discussion

In this section we give the CFDRM and Keller-box method results for the two examples described above. The Keller-box method is an implicit finite difference scheme which is made up of four major steps. First, the governing differential equation is reduced to a system of first-order equations. The resulting first-order equations are discretized using central finite differences. The resulting algebraic equations are linearised by Newton's method and written in matrix vector form. The linear system obtained is then solved by block tridiagonal elimination technique.

Example 1 was solved for M=1 . The value of _η∞ , which correspond to η=∞ , was approximated to be equal to 10 . This value was found to be large enough to give consistent results. For the CFDRM, uniform stepsizes Δη=0.1 in the η -direction and Δξ=0.0025 in the ξ -direction were used to get results accurate to six decimal places. For the Keller-box, stepsizes Δη=0.005 and Δξ=0.004 were used. Table 1 shows the skin friction ^{f[variant prime][variant prime]} (0,ξ) for Example 1 obtained using both CFDRM and Keller-box method. The times taken to compute the solution at different values of ξ are also given in Table 1. It can be seen that the CFDRM takes a very short period of time to compute the solution compared to the Keller-box method. The Keller-box required quite a small Δη to get a solution accurate to six decimal places. This means that more points are needed in the Keller-box compared to the CFDRM. Figure 1 shows the effects of varying the magnetic field parameter (M ) on the skin friction ^{f[variant prime][variant prime]} (0,ξ) . The parameters M=1 , c=0.5 , and Pr=0.7 were used to solve the governing equations in Example 2. The value of _η∞ was approximated to be equal to 10 . For the CFDRM, uniform stepsizes Δη=0.1 and Δξ=0.0025 to get results accurate to six decimal places. For the Keller-box, stepsizes Δη=0.005 and Δξ=0.002 were used. The quantities of physical interest (^{f[variant prime][variant prime]} (0,ξ) , ^{s[variant prime][variant prime]} (0,ξ) , and ^{g[variant prime]} (0,ξ) ) are tabulated for both CFDRM and Keller-box and are shown in Table 2. The results are in good agreement with six decimal places as depicted in the table. Table 3 shows the times taken by each method at various values of ξ . It is evident from both Tables 2 and 3 that the CFDRM computes the solution in a very short period of time compared to the Keller-box. The main reason is that a very small Δη and hence a large number of points are required in the Keller-box method for the required accuracy. Graphical solutions of ^{f[variant prime][variant prime]} (0,ξ) , ^{s[variant prime][variant prime]} (0,ξ) , and ^{g[variant prime]} (0,ξ) in the entire range 0...4;ξ...4;1 for Example 2 are shown in Figures 2, 3, and 4. The solution is plotted for different values of M and c . Comparing the plots with those of [8], it can be seen that the results of the CFDRM are comparable with their homotopy analysis method results. Xu et al. [8] used the homotopy analysis method to get series solutions which are accurate and uniformly valid for all dimensionless time 0...4;ξ...4;1 . From these figures it is clear that the CFDRM is able to give results which are accurate and uniformly valid for all dimensionless time 0...4;ξ...4;1 .

Table 1: CFDRM and Keller-box results for the skin friction ^{f[variant prime][variant prime]} (0,ξ) in Example 1 when M=1 .

ξ	CFDRM	Keller-box
f [variant prime] [variant prime] ( 0 , ξ )	Time (sec)	f [variant prime] [variant prime] ( 0 , ξ )	Time (sec)
0.1	-0.665032	0.234	-0.665031	3.427
0.3	-0.854591	0.687	-0.854591	10.029
0.5	-1.029610	1.156	-1.029610	16.657
0.7	-1.191787	1.606	-1.191787	23.267
0.9	-1.342594	2.110	-1.342594	31.140

Table 2: CFDRM and Keller-box results for ^{f[variant prime][variant prime]} (0,ξ),^{s[variant prime][variant prime]} (0,ξ) , and ^{g[variant prime]} (0,ξ) in Example 2 when M=1 , c=0.5 , and Pr=0.7 .

ξ	f [variant prime] [variant prime] ( 0 , ξ )	s [variant prime] [variant prime] ( 0 , ξ )	g [variant prime] ( 0 , ξ )
CFDRM	Keller-box	CFDRM	Keller-box	CFDRM	Keller-box
0.1	-0.674444	-0.674444	-0.327586	-0.327586	-0.483085	-0.483085
0.3	-0.880408	-0.880408	-0.414565	-0.414565	-0.503506	-0.503506
0.5	-1.068930	-1.068930	-0.496369	-0.496369	-0.521021	-0.521021
0.7	-1.242068	-1.242068	-0.573224	-0.573224	-0.534064	-0.534064
0.9	-1.401598	-1.401598	-0.645382	-0.645382	-0.537130	-0.537130
Time (sec)	6.074	238.769

Table 3: Computational times for Example 2 at various values of ξ .

ξ	CFDRM	Keller-box
0.1	0.635	22.644
0.3	1.802	65.809
0.5	3.032	109.886
0.7	4.256	157.737
0.9	5.559	215.769

Figure 1: CFDRM solution of -^{f[variant prime][variant prime]} (0,ξ) for different values of M for Example 1.

[figure omitted; refer to PDF]

Figure 2: CFDRM solution of -^{f[variant prime][variant prime]} (0,ξ) and -^{s[variant prime][variant prime]} (0,ξ) for different values of M when c=0.5 and Pr=0.7 .

[figure omitted; refer to PDF]

Figure 3: CFDRM solution of -^{f[variant prime][variant prime]} (0,ξ) and -^{s[variant prime][variant prime]} (0,ξ) for different values of c when M=1 and Pr=0.7 .

[figure omitted; refer to PDF]

Figure 4: CFDRM solution of -^{g[variant prime]} (0,ξ) for different values of M and c , respectively, when Pr=0.7 .

[figure omitted; refer to PDF]

6. Conclusion

In this paper, two unsteady boundary layer flow problems have been investigated. The nonlinear partial differential equations governing the unsteady boundary layer flow of the two examples considered are solved using the compact finite difference relaxation method (CFDRM). The CFDRM had previously been applied only to ordinary differential equations modelling boundary layer problems [17]. This work has successfully extended its application to partial differential equations. The results obtained are in good agreement with results obtained using the Keller-box. The CFDRM is computationally faster than the Keller-box method. This is because the CFDRM gives highly accurate solutions on coarser grids.

References

[1] B. C. Sakiadis, "Boundary-layer behavior on continuous solid surfaces: II. The boundary layer on a continuous flat surface," AIChE Journal , vol. 7, no. 2, pp. 221-225, 1961.

[2] J. C. Williams III, T. B. Rhyne, "Boundary layer development on a wedge impulsively set into motion," SIAM Journal on Applied Mathematics , vol. 38, no. 2, pp. 215-224, 1980.

[3] T. Cebeci, P. Bradshaw Physical and Computational Aspects of Convective Heat Transfer , pp. xii+487, Springer, New York, NY, USA, 1984.

[4] R. Seshadri, N. Sreeshylan, G. Nath, "Unsteady mixed convection flow in the stagnation region of a heated vertical plate due to impulsive motion," International Journal of Heat and Mass Transfer , vol. 45, no. 6, pp. 1345-1352, 2002.

[5] R. Nazar, N. Amin, I. Pop, "Unsteady boundary layer flow due to a stretching surface in a rotating fluid," Mechanics Research Communications , vol. 31, no. 1, pp. 121-128, 2004.

[6] R. Nazar, N. Amin, D. Filip, I. Pop, "Unsteady boundary layer flow in the region of the stagnation point on a stretching sheet," International Journal of Engineering Science , vol. 42, no. 11-12, pp. 1241-1253, 2004.

[7] S. Liao, "An analytic solution of unsteady boundary-layer flows caused by an impulsively stretching plate," Communications in Nonlinear Science and Numerical Simulation , vol. 11, no. 3, pp. 326-339, 2006.

[8] H. Xu, S.-J. Liao, I. Pop, "Series solutions of unsteady three-dimensional MHD flow and heat transfer in the boundary layer over an impulsively stretching plate," European Journal of Mechanics B , vol. 26, no. 1, pp. 15-27, 2007.

[9] X. You, H. Xu, I. Pop, "Homotopy analysis of unsteady heat transfer started impulsively from rest along a symmetric wedge," International Communications in Heat and Mass Transfer , vol. 37, no. 1, pp. 47-51, 2010.

[10] T. Fan, H. Xu, I. Pop, "Unsteady stagnation flow and heat transfer towards a shrinking sheet," International Communications in Heat and Mass Transfer , vol. 37, no. 10, pp. 1440-1446, 2010.

[11] Z. Abbas, T. Hayat, M. Sajid, S. Asghar, "Unsteady flow of a second grade fluid film over an unsteady stretching sheet," Mathematical and Computer Modelling , vol. 48, no. 3-4, pp. 518-526, 2008.

[12] I. Ahmad, M. Sajid, T. Hayat, M. Ayub, "Unsteady axisymmetric flow of a second-grade fluid over a radially stretching sheet," Computers & Mathematics with Applications , vol. 56, no. 5, pp. 1351-1357, 2008.

[13] T. Hayat, Z. Abbas, "Heat transfer analysis on the MHD flow of a second grade fluid in a channel with porous medium," Chaos, Solitons & Fractals , vol. 38, no. 2, pp. 556-567, 2008.

[14] A. Ali, A. Mehmood, "Homotopy analysis of unsteady boundary layer flow adjacent to permeable stretching surface in a porous medium," Communications in Nonlinear Science and Numerical Simulation , vol. 13, no. 2, pp. 340-349, 2008.

[15] A. Mehmood, A. Ali, T. Shah, "Heat transfer analysis of unsteady boundary layer flow by homotopy analysis method," Communications in Nonlinear Science and Numerical Simulation , vol. 13, no. 5, pp. 902-912, 2008.

[16] A. Alizadeh-Pahlavan, K. Sadeghy, "On the use of homotopy analysis method for solving unsteady MHD flow of Maxwellian fluids above impulsively stretching sheets," Communications in Nonlinear Science and Numerical Simulation , vol. 14, no. 4, pp. 1355-1365, 2009.

[17] P. G. Dlamini, S. S. Motsa, M. Khumalo, "On the comparison between compact finite difference and pseudospectral approaches for solving similarity boundary layer problems," Mathematical Problems in Engineering , vol. 2013, 2013.

[18] S. K. Lele, "Compact finite difference schemes with spectral-like resolution," Journal of Computational Physics , vol. 103, no. 1, pp. 16-42, 1992.

[19] M. Sari, G. Gürarslan, "A sixth-order compact finite difference scheme to the numerical solutions of Burgers' equation," Applied Mathematics and Computation , vol. 208, no. 2, pp. 475-483, 2009.

[20] P.-G. Zhang, J.-P. Wang, "A predictor-corrector compact finite difference scheme for Burgers' equation," Applied Mathematics and Computation , vol. 219, no. 3, pp. 892-898, 2012.

[21] A. Shah, L. Yuan, A. Khan, "Upwind compact finite difference scheme for time-accurate solution of the incompressible Navier-Stokes equations," Applied Mathematics and Computation , vol. 215, no. 9, pp. 3201-3213, 2010.

[22] J. Li, M. R. Visbal, "High-order compact schemes for nonlinear dispersive waves," Journal of Scientific Computing , vol. 26, no. 1, pp. 1-23, 2006.

[23] B. Düring, M. Fournié, A. Jüngel, "High order compact finite difference schemes for a nonlinear Black-Scholes equation," International Journal of Theoretical and Applied Finance , vol. 6, no. 7, pp. 767-789, 2003.

[24] M. Bastani, D. K. Salkuyeh, "A highly accurate method to solve Fisher's equation," Pramana , vol. 78, no. 3, pp. 335-346, 2012.

[25] M. Sari, "Solution of the porous media equation by a compact finite difference method," Mathematical Problems in Engineering , vol. 2009, 2009.

[26] M. Sari, G. Gürarslan, "A sixth-order compact finite difference method for the one-dimensional sine-Gordon equation," International Journal for Numerical Methods in Biomedical Engineering , vol. 27, no. 7, pp. 1126-1138, 2011.

[27] A. H. Srinivasa, A. T. Eswara, "Unsteady MHD laminar boundary layer flow due to an impulsively stretching surface," in Proceedings of the World Congress on Engineering (WCE '11), vol. 1, pp. 252-255, London, UK, July 2011.

[28] S. S. Motsa, Z. G. Makukula, "On spectral relaxation method approach for steady von Kármán flow of a Reiner-Rivlin fluid with Joule heating, viscous dissipation and suction/injection," Central European Journal of Physics , vol. 11, no. 3, pp. 363-374, 2013.

[29] S. S. Motsa, P. Dlamini, M. Khumalo, "A new multistage spectral relaxation method for solving chaotic initial value systems," Nonlinear Dynamics , vol. 72, no. 1-2, pp. 265-283, 2013.

[30] S. S. Motsa, P. G. Dlamini, M. Khumalo, "Solving hyperchaotic systems using the spectral relaxation method," Abstract and Applied Analysis , vol. 2012, 2012.