Abstract

This paper introduces two novel numerical algorithms for the efficient solution of coupled systems of nonlinear boundary value problems. The methods are benchmarked against existing methods by finding dual solutions of the highly nonlinear system of equations that model the flow of a viscoelastic liquid of Oldroyd-B type in a channel of infinite extent. The methods discussed here are the spectral relaxation method and spectral quasi-linearisation method. To verify the accuracy and efficiency of the proposed methods a comparative evaluation of the performance of the methods against established numerical techniques is given.

1. Introduction

Exact solutions to a wide class of problems in engineering and science are generally available only for a limited range of problems. For this reason the quest for new techniques and the improvement of existing techniques for finding solutions of nonlinear equations are an ongoing challenge in engineering and science. In addition to the classical numerical methods, such as those based on finite differences, finite elements, and finite volume techniques, there is currently a wide variety of methods for nonlinear equations, such as, among others, linearisation methods [14] and the transform methods of Fokas [57]. Advances in decomposition and variational methods in recent years have added to the repertoire of available techniques for finding solutions of nonlinear BVPs. In recent years these methods have been further augmented by He’s [8, 9] and Liao’s [10, 11] homotopy based methods and their various variants such as the spectral-homotopy analysis of Motsa et al. [12, 13] and Liao’s [14] generalized boundary element method.

This paper introduces two novel techniques based on a combination of linearisation techniques and spectral methods and that allow for simple and straightforward integration of systems of ordinary differential equations on finite and infinite domains. We present an overview of these techniques and provide a comparative evaluation of the two methods against results in the literature.

A recent study by Bachok et al. [15] investigated the two-component convection in a viscoelastic liquid of Oldroyd-B type occupying a horizontal channel of infinite extent and depth . The flow is governed by the coupled system of equations where , , and represent, respectively, the nondimensional stream function, temperature, and nanoparticle volume fraction, is the Prandtl number, is the Lewis number, is the Brownian motion number, and and are constant dimensionless parameters.

Equations (1)–(3) are solved subject to the boundary conditions where is the stretching/shrinking velocity and is a constant mass flux. In Bachok et al. [15] it was shown that dual solutions of (1)–(3) exist for , where    is some critical value of . No solutions exist for .

In the present study we revisit the solution of the system of nonlinear equations (1)–(5) using the spectral relaxation method (Motsa [16]) and the spectral quasi-linearization method (SQLM). The objective of this study is to give a comparative analysis of the performance of the two techniques in finding solutions of coupled highly nonlinear problems in fluid mechanics. We determine, inter alia, the accuracy of each method and demonstrate how dual solutions can be obtained using the SQLM. In addition, we present a systematic approach of obtaining critical parameter values and multiple solutions of the governing equations.

2. The Spectral Relaxation Method

In this section we describe the development of the spectral relaxation method (SRM) for the solution of the nonlinear system (1)–(3). The system of equations is decoupled and the resulting subsystems are solved in a sequential manner. The method appears to be particularly effective for nonlinear systems of differential equations in which some of the unknown functions have exponentially decaying profiles. The algorithm for the method when applied to the system (1)–(3) may be summarized as follows.(1)Reduce the order of the momentum equation from three to one by using the transformation (2)Assume that is known from some initial guess and arrange the transformed equations in a particular order, placing the equations with the least unknowns at the top of the equations list. This gives with (2) and (3) placed below (8) in their original form.(3)Assign the labels , and to the order list of (7), (8), (2), and (3), respectively.(4)In the equation for (1st equation), the iteration scheme is developed by assuming that only linear terms in are to be evaluated at the current iteration level (denoted by ) and all other terms (linear and nonlinear) in , and are assumed to be known from the previous iteration (denoted by ). In addition nonlinear terms in are also evaluated at the previous iteration.(5)In developing the iteration scheme for the next equation , only linear terms in are evaluated at the current iteration level with all other terms evaluated at the previous level, except which is now known from the solution of the first equation.(6)This process is repeated in the equations for using the updated solutions for . The same procedure is effected on the equation for , now using the updated solutions for , and .

The strategy used for decoupling the system of equations is analogous to the Gauss-Seidel relaxation method which is normally used in solving linear algebraic system of equations. In the context of the algorithm described above, we obtain the following iteration scheme: subject to the boundary conditions Given a set of initial approximations , (9) can be solved in a sequential manner, one after the other, using an appropriate numerical method. It must be noted that (9) now forms a sequence of linear differential equations with variable coefficients which can easily be solved using spectral collocation methods. The spectral methods are preferred here because of their remarkably high accuracy and relative simplicity in discretizing and subsequent solution of variable-coefficient linear differential equations with smooth solutions over simple domains. Spectral methods, such as the Chebyshev spectral collocation methods, have been found to be very efficient in discretizing and solving other iteration schemes for solving nonlinear equations of boundary layer type. Examples of recently developed spectral based method for solving boundary layer type equations include the spectral homotopy analysis method [12, 13] and the successive linearisation method [17, 18]. The SLM has been successfully applied to different nonlinear boundary value problems in fluid mechanics; see Awad et al. [19], Motsa and Sibanda [20], and Sibanda et al. [17, 18]. A large volume of literature exists on the practical implementation of spectral collocation methods including the books [21, 22] which the interested reader may find useful.

In applying the spectral collocation method, we find the unknown function at collocation points by requiring that (9) be satisfied exactly at these points. A convenient set of collocation points is the Gauss-Lobatto points defined on by We approximate the derivatives of the unknown functions, say , using the so-called differentiation matrix of size which is computed as at the collocation points as the matrix vector product where is a finite value used to numerically approximate the conditions at infinity, , and is the vector function at the collocation points. The variable is defined from the linear transformation used to convert the truncated interval into the interval on which the spectral method can be implemented. The second order derivatives are obtained as powers of ; that is Applying the Chebyshev spectral collocation method on (9) we obtain the following sequence of matrix equations: where where is an diagonal matrix, diag denotes a diagonal matrix, and , , , and correspond to the approximation of the functions , , , and when evaluated at the collocation points. The SRM iteration scheme (14)-(15) is solved for , starting from suitable initial guesses , and . We choose the following initial approximations which satisfy the boundary conditions of the governing equations:

3. The Spectral Quasi-Linearisation Method

In this section, we present an alternative method for solving the governing system of (1)–(3). We employ the quasi-linearization method (QLM) which is a generalization of Newton-Raphson method. This method was initially proposed by Bellman and Kalaba [1] for solving nonlinear boundary value problems. Here we extend the application of the QLM to the governing equations (1)–(3) and discretize the QLM equations using the spectral collocation method as described in the previous section. In this work, the spectral method based quasi-linearization method is referred to as the spectral quasi-linearisation method (SQLM).

To develop the SQLM scheme it is convenient to write the governing system as a sum of its linear and nonlinear components as for , where (1), (2), and (3) are labeled using , and 3, respectively. Thus, we have

The QLM is equivalent to multivariable Taylor series expansion of the nonlinear functions assuming that the difference in the values of each unknown function at the current iteration denoted by, say, and the previous iteration is small. Thus, the QLM scheme corresponding to (1)–(3) is given by subject to

To derive the above QLM scheme, it was noted that the solution for can be resolved independently of (2) and (3). Thus, in deriving the iteration schemes for and , the solution for at the current iteration () was assumed to be known.

Applying the Chebyshev pseudo-spectral method on (22) we obtain the following SQLM scheme in matrix form: subject to the boundary conditions where where denotes an matrix of zeros and denotes diag( ); that is, the vector elements are placed on the main diagonal of a matrix whose entries everywhere else are zero. Starting from the initial approximations (19), the approximate SQLM solutions for , and are obtained by solving (24) subject to the boundary conditions (25).

4. Determining System Critical Values

The system (1)–(3) was reported in Bachok et al. [15] to have multiple solutions and bifurcation curves on the planes , , and . The parameters , , and are, respectively, skin friction coefficient, local Nusselt number, and the local Sherwood number and are defined in terms of the governing parameters as

In this section, we demonstrate how the spectral quasi-linearisation method (SQLM) may be used to determine the critical values of beyond which solutions do not exist.

To determine the critical values , we regard and to be additional unknown variables. We define and differentiate (1) with respect to to obtain subject to where . Criticality of the curve corresponds to ; thus, . Applying quasi-linearisation on (1) and (28)-(29) gives subject to Equations (30) and (31) form a system of four equations which are solved iteratively for the unknowns starting from suitable initial approximations . Applying the spectral collocation method on (30) and (31) gives subject to the boundary conditions where In matrix form, (32) and (33) can be expressed as

5. Results and Discussion

The governing systems of (1) and (3) were solved using the spectral relaxation method (SRM) and spectral quasi-linearisation method (SQLM). In this section, we present the results of the numerical computations for the velocity, temperature, and concentration profiles for the various input parameters. The convergence and stability of the iteration schemes can be assessed by considering the norm of the difference in the values of functions between two successive iterations. For a converging iteration scheme, the difference in the norm of the value of the functions computed at successive iterations is expected to decrease with an increase in the number of iterations. The accuracy of the present results was verified by comparing with the results of Bachok et. al. [15]. In the calculations presented here, the values of the governing physical parameters were chosen deliberately to match the results of [15] in order to enable effective comparison. The effect of the number of collocation points was examined in order to select the smallest value of which gives a consistent solution. This is achieved by repeatedly solving the governing equations using the proposed iteration schemes with different values of until the consistent solution within a small error tolerance level was reached. In this work, the error tolerance was set to be . Unless otherwise specified, the number of collocation points used in this study is . The , that is, at infinity, was also chosen in such a way that further changes in its magnitude do not produce changes in the values of and its derivatives within the tolerance level . A value of was found to be appropriate in all cases discussed in this investigation.

Figures 1 and 2 show the variation of the error , the difference of the norm of computed using the spectral relaxation method at different iterations for different system parameter values. In all three cases it can be seen that decreases as the number of iterations increases. This shows that the SRM converges for the selected choice of parameters. The convergence of the SRM improves with an increase in the unsteadiness parameter and mass suction parameter but decreases with an increase in the stretching parameter .

Figure 3 gives the variation of the SQLM generated error with the number of iterations. The convergence of the SQLM is seen to be rapid for the first few iterations then plateaus after about five iterations. This result presents the main difference between the SRM and SQLM. The SRM results in Figures 1 and 2 indicate strict convergence of the method with the error progressively reduced with an increase in the number of iterations. For the same number of collocations points, the SQLM will converge faster than the SRM but is less accurate. This observation suggests that the convergence rate of the SQLM approach is superior. However, it turns out that the SRM yields a better accuracy for the three ODEs system.

To obtain the multiple solutions we choose an initial guess of the form where is an unknown constant which must be carefully selected. It is observed that when varying between negative and positive values, the concavity of the velocity profile changes in the region near . By fixing the values of the governing physical parameters and varying the values of in the SQLM implementation, the values of that give the multiple solutions can be identified from the maximum of the residual of the momentum equation (1). Figure 4 shows the variation of the maximum residual of the momentum equation (1) against at different SQLM iterations when . It can be seen that the graph develops two local minima. The multiple solutions were obtained by setting to be equal to the value of at which the local minima occurs. For example, when , two solutions were obtained when and . It is worth mentioning here that the SRM approach did not yield the second solution.

Figure 5 shows the two solution profiles and , respectively, for different values of and fixed values of and . The corresponding temperature and concentration profiles are depicted in Figure 6.

In Table 1 we give a comparison between the present SQLM results and the results of [15] for the critical values of for different values of and . We note that the present results are in good agreement with the results reported in Bachok et al. [15]. This shows the accuracy of the SQLM when compared to the shooting method as used in [15].

Figure 7 shows the variation of the local skin friction and local Nusselt and Sherwood numbers with for different values of and and corresponding to the other parameters values reported in Bachok et al. [15]. The SQLM easily resolves the dual solutions when . These results obtained using the SQLM are qualitatively and quantitatively the same as those of Bachok et al. [15]. The results show the two solutions branches for each value of the unsteadiness parameter and the mass flux parameter .

6. Conclusion

In this study we have used two methods, the spectral relaxation method and the spectral quasi-linearisation method to solve the highly nonlinear equations that describe the unsteady heat transfer in a nanofluid over a permeable stretching or shrinking surface. The results obtained by Bachok et al. [15] have been confirmed. In terms of the accuracy, convergence rates, and robustness of the techniques used, we have demonstrated that(i)The SQLM converges faster than SRM;(ii)The SRM is more accurate than the SQLM;(iii)The SQLM can resolve multiple solutions and SRM only yields one solution.

Conflict of Interests

The authors declare that there is no conflict of interests.

Acknowledgments

The work of S. S. Motsa is based on the research supported in part by the National Research Foundation of South Africa (Grant no. 85596).