The Effect of Thermophoresis on Unsteady Oldroyd-B Nanofluid Flow over Stretching Surface

There are currently only a few theoretical studies on convective heat transfer in polymer nanocomposites. In this paper, the unsteady incompressible flow of a polymer nanocomposite represented by an Oldroyd-B nanofluid along a stretching sheet is investigated. Recent studies have assumed that the nanoparticle fraction can be actively controlled on the boundary, similar to the temperature. However, in practice, such control presents significant challenges and in this study the nanoparticle flux at the boundary surface is assumed to be zero. We have used a relatively novel numerical scheme; the spectral relaxation method to solve the momentum, heat and mass transport equations. The accuracy of the solutions has been determined by benchmarking the results against the quasilinearisation method. We have conducted a parametric study to determine the influence of the fluid parameters on the heat and mass transfer coefficients.


Introduction
A wide variety of fluids, such as polymer solutions, plastics, pulps, emulsions, blood plasma, chocolate, tomato sauce, mustard, mayonnaise, toothpaste, asphalt, some greases and sewage, petroleum products, oils, etc., are non-Newtonian in character. Such fluids have a non-linear shear stress-strain rate relationship. The equations that model the flow of these fluids are generally of a higher order than the Navier-Stokes equations. During the last several decades, the study of flow of non-Newtonian fluids, has received considerable amount of research interest, due to the relevance of non-Newtonian flows to a large number of engineering and manufactural applications such as in the processing of synthetic fibers, food, polymers melts and pharmaceutical products.
The unsteady boundary layer flow of non-Newtonian fluids due to stretching surface is an important field of research in fluid mechanics. In many practical problems, the flow could be unsteady due to a time dependent free stream velocity. There are several transport processes with surface mass transfer where the buoyancy force arises from thermal diffusion caused by the temperature gradient. It is interesting as well as useful to investigate the combined effects of thermal diffusion and surface mass transfer on a stretching surface where the free stream velocity varies arbitrarily with time. The stretched boundary layer flow with heat transfer has numerous applications in engineering and industrial processes such as when polymer sheets are extruded continuously from a dye, in the annealing and thinning of copper wires, paper production and glass blowing, aerodynamic extrusion of plastic sheets and the cooling of infinite metallic plates in cooling baths. Mixed convection or buoyancy driven flows over a surface occur naturally in geothermal and petroleum recovery processes, solid matrix heat exchanges, the cooling of nuclear reactors and other practical problems. Historically, the study of a boundary layer flow along a stretching surface has its origins in the work of Sakiadis [1]. Tsou et al. [2] presented a combined experimental and analytical study of the stretching flow, which, in essence demonstrated that such a flow is physically possible. Crane [3] further extended the Sakiadis study to a linearly stretching plate in a quiescent fluid and presented an exact analytic solution. Subsequent work has looked at various aspects of the stretching sheet problem, but most have been concerned with how the flow is affected by or responds to changes in various fluid and surface parameters, see for instance, [4][5][6][7][8][9].
Conventional fluids such as water, ethylene glycol and oil have low heat transfer characteristics owing to their low thermal conductivity. A recent technique to enhance the thermal conductivity of these base fluids is to suspend nano-sized metallic particles such as aluminum, titanium, gold, copper, iron or their oxides in the conventional base fluids resulting in what has come to be name as a "nanofluid," Choi [10]. Nanofluids have considerably improved thermophysical properties such as thermal conductivity, thermal diffusivity, viscosity and convective heat transfer coefficient compared to base fluids. In the last few years the flow of nanofluid through different geometries, and under various physical assumptions have been studied by several authors such as [11][12][13][14].
The Oldroyd-B constitutive model describes a subclass of non-Newtonian fluids that adequately describe the behaviour of some viscoelastic fluids such as dilute polymer solutions. The Oldroyd-B fluid can describe stress-relaxation, creep and the normal stress differences but it cannot describe either shear thinning or shear thickening, a phenomenon that is exhibited by many polymer materials. Nonetheless, this model is perhaps one of the most successful models for describing the response of some polymeric liquids [15]. Some investigations of Oldroyd-B fluids have been done, by, among others, Hayat et al. [16] who presented a study of the threedimensional flow of an Oldroyd-B fluid due to a stretching surface with convective boundary conditions. Siddiqui et al. [17] investigated the unsteady flow of an incompressible Oldroyd-B fluid between two infinite parallel plates subject to slip between the plates and the fluid. Jamil et al. [18] further studied the unsteady flow of an Oldroyd-B fluid and solved the model equations using finite Hankel transforms. Sohail et al. [19] investigated the two-dimensional steady incompressible Oldroyd-B nanofluid flow past a stretching sheet. The thermophoresis and radiation effects on Heat and mass transfer characteristics in three-dimensional flow of an Oldroyd-B fluid due to a bi-directional stretching surface were investigated by Shehzad et al. [20]. Related studies include, among others, [21][22][23][24][25]. Excellent survey papers on Oldroyd-B fluids can be found in studies by Rajagopal and his colleagues [26][27][28].
There are currently only a few theoretical studies on the convective boundary layer flow and heat transfer in an Oldroyd-B nanofluid over a stretching surface. Including nanoparticles in dilute polymer solutions imparts improved thermophysical properties to the polymer materials. This could, for example, include improved electrical, mechanical and optical properties, Chao [29]. Recent exception is the study by Khan et al. [30]. These studies all assumed that the nanoparticle volume fraction at the boundary is actively controlled. To the best of the authors' knowledge, there are as yet no studies of convective heat transport in polymer nanocomposites such as an Oldroyd-B nanofluid flow in which the nanofluid particle fraction on the boundary is not actively managed. The objective of this study is therefore to investigate thermophoresis effects in unsteady Oldroyd-B nanofluid flow along a stretching surface. The problem is formulated under the assumption that the nanoparticle volume fraction at the boundary is not actively controlled, see Nield and Kuznetsov [31][32][33][34]. The highly non-linear momentum, heat and mass transfer equations are solved numerically using the spectral relaxation and quasi-linearization methods, see Motsa et al. [35,36].

Mathematical Formulation
Consider the unsteady two-dimensional Oldroyd-B nanofluid flow over a stretching surface. The sheet stretches along the plane y = 0 with the flow confined in the region y > 0. The nanoparticle flux at the boundary surface is assumed to be zero. The surface is stretched with a velocity U = bx/(1−at) where a and b are positive constants. Both the nanofluid and the surface are kept at a constant temperature T w where T w > T 1 is for a heated surface and T w < T 1 corresponds to a cooled surface. The geometry of the problem is shown in Fig 1. Applying the Oberbeck-Boussinesq and the boundary layer approximations to the basic equations of an incompressible non-Newtonian fluid, we obtain The Cauchy stress tensor τ and extra stress tensor S are defined as where D Dt ¼ @ðÁÞ @t þ ðV Á rÞðÁÞ À ðrVÞðÁÞ À ðÁÞðrVÞ Ã ; is the covariant differentiation, the Ã denotes the matrix transpose,d 1 ðtÞ andd 2 ðtÞ are the relaxation and retardation times respectively. The velocity field is V ¼ ðuðx; y; tÞ; vðx; y; tÞÞ; where u and v are the velocity components along the x and y directions respectively. The first Rivlin-Ericksen tensor A is defined as The governing equations take the form (see Nadeem et al. [19]) subject to the boundary conditions where ν is the kinematic viscosity, T and are the local fluid temperature and concentration volume fraction, T 1 and 1 are the fluid temperature and ambient concentration volume fraction respectively, α m is the effective thermal diffusivity, D B is the Brownian diffusion coefficient, D T is thermophoresis diffusion coefficient, τ = (ρc) f /(ρc) P is the ratio between the effective heat capacity of the nanoparticle material and heat capacity of the fluid. We introduce the following similarity transformations cðx; y; tÞ ¼ where θ and ϕ are the dimensionless temperature and nanoparticle volume fraction respectively and η is a similarity variable. The physical stream function ψ(x, y, t) automatically ensures that mass conservation given in Eq (6). The velocity components are readily obtained as where f 0 is the dimensionless velocity. Eqs (7)- (8) can be presented as with the boundary conditions where S is the dimensionless measure of the unsteadiness, β 1 and β 2 are the Deborah numbers in terms of relaxation and retardation times, respectively, the Prandtl number Pr, the Brownian motion parameter N b , the thermophoresis parameter N t , the Lewis number Le. These parameters are defined as The non-dimensional form of the Nusselt number and Sherwood number that describe the heat and nanoparticle volume fraction transfer at the surface are Nu x =Re 1 2 x ¼ Ày 0 ð0Þ and Shw x =Re where Re x = Ux/ν is the local Reynolds number.

Method of Solution
To solve Eqs (12)-(14) along with the boundary conditions Eqs (15)-(16), the spectral relaxation method (SRM) was used, see Motsa et al. [37][38][39]. This method is preferred since it has been shown to be accurate and generally easier to use compared to other common numerical methods such as finite differences.
We start by reducing the order of Eq (12) from fourth to third order. To this end, we set f 0 = g, so that Eq (12) becomes The spectral relaxation algorithm decouples the system of governing Eqs (12)- (14). From the decoupled equations an iteration scheme is developed by evaluating linear terms at the current iteration level r + 1 and the nonlinear terms at the previous iteration level r. Applying the SRM to Eqs (13)- (14) and (17)-(17) gives the following linear ordinary differential equations; where a i , b i and c i (i = 1, 2, . . .) are given by Starting from given initial approximations f 0 , g 0 , θ 0 and ϕ 0 , Eqs (18)-(21) can be solved iteratively using any suitable numerical method. We opt to use the spectral collocation methods for its accuracy. We find the unknown function at collocation points by requiring that Eqs (18)- (21) be satisfied exactly at these points. A convenient set of collocation points is the Gauss-Lobatto points defined by For convenience, in numerical computations the semi-infinite domain is approximated by the truncated domain [0, L]. Then using the linear transformation η = L(ω + 1)/2, we convert [0, L] into the interval [−1, 1] in which the spectral method can be used, where L = η 1 is a finite number selected to be large enough to represent the behaviour of the flow properties when η is very large. The derivatives are defined as where N + 1 is the number of collocation points, D = 2D/L and f = [f(ω 0 ), f(ω 1 ), . . ., f(ω N )] T is the vector of unknown functions at the collocation points. Applying the Chebyshev spectral collocation method to the system Eqs (18)-(21), we obtain the following matrix equations Here I is an (N + 1) × (N + 1) diagonal matrix, diag[Á] denotes a diagonal matrix. We choose suitable initial guesses f 0 , g 0 , θ 0 and ϕ 0 which satisfy the boundary conditions of governing equations as

Results and Discussion
Eqs (12)- (14) along with the boundary conditions Eqs (15)-(16), were solved numerically using both the spectral relaxation method (SRM) and the quasi-linearization method (QLM), see Bellman and Kalaba [40]. Here the QLM has been used as a benchmarking tool to test the accuracy, and hence the reliability of the SRM results. The results showing the effects of various parameters on the skin-friction coefficient and the heat transfer rate on the unsteady Oldroyd-B nanofluid are given in Tables 1-3.  Table 1 gives a comparison between the present results and the results obtained by Sharidan [41] and Pal [42] for the skin friction. There is a good agreement between the two sets of results with the SRM having converged at the fifth order up to six decimal places. The SRM results are further validated by comparison with results generated using a quasilinearisation method. The quasilinearisation method has also been used recently, albeit in a slightly different form, by Ibrahim and Shanker [43] Table 2 shows a further comparison of the spectral relaxation and quasilinearisation results for the skin friction coefficient while Table 3 shows the variation of the heat transfer rate for different values of dimensionless unsteadiness parameter and Deborah numbers. The comparison of the two methods shows an excellent agreement between the numerical results obtained by the spectral relaxation and the quasilinearisation methods. In addition, Table 2 shows that the skin friction coefficient increases with increasing values of the unsteadiness parameter S and the relaxation time in terms of the Deborah number β 1 and decreases with increasing retardation time or Deborah number β 2 . The heat transfer coefficient −θ 0 (0) is however shown to decrease with β 1 in Table 3. he increase in the skin friction coefficient with the flow unsteadiness has also been observed in earlier studies, such as in Ibrahim and Shanker [43] and Mukhopadhyay et al. [44]. It has been suggested in Mukhopadhyay et al. [44] that a decrease in the skin friction coefficient may be important in coating processes where higher stretching speeds may be achieved for smaller pulling forces. The study by Mukhopadhyay et al. [44] was on unsteady flow in a Casson fluid and further showed that for Casson fluids, the temperature decreased significantly with unsteadiness.
We observe further that heat transfer coefficient increases with increased unsteadiness and β 2 . The negative values of the nanoparticle profile are due to the fact that the effect of thermophoresis is such that an elevation above the ambient surface temperature leads to a reduction in the relative value of the nanoparticle fraction at the surface (see Kuznetsov and Nield [34]). Figs 2 and 3 show the velocity profiles for different values of the unsteadiness parameter. We observe that the velocity distribution and the momentum boundary layer thicknesses reduce with an increase in the unsteadiness parameter. This finding is in line with the earlier findings of Ibrahim and Shanker [43], and shows that even in the absence of an applied magnetic field, the velocity profiles decrease with the unsteadiness parameter. Figs 4 and 5 show the dimensionless temperature and concentration volume fraction profiles respectively for selected values of S. The steepness in both the temperature and concentration profiles decreases reducing the thicknesses of both the thermal and concentration volume fraction boundary layers. These results also follow a similar trend as observed by Ibrahim and Shanker [43] and Mukhopadhyay et al. [44]. Further, it may be pointed out that the concentration volume fraction increases from negative to positive values until boundary layer separation occurs. The magnitude of the volume fraction concentration increases up to a critical point and then decreases to zero. Figs 6 and 7 show the influence of the Deborah numbers β 1 and β 2 on the velocity profiles. Higher Deborah numbers are indicative that the Oldroyd-B nanofluid is stretched. The nanofluid velocity f 0 (η) and the momentum boundary layer thickness decrease with increasing β 1 which is not an unexpected result since it is well known that the viscoelastic fluid resists the motion of the fluid. Figs 8 and 9 show the effects of β 1 on the temperature and concentration profiles, respectively. As β 1 increases, both the nanofluid temperature and the concentration volume fraction increase enhancing both the thermal and the concentration boundary layer thicknesses. Table 3. Comparison of heat transfer rate −θ 0 (0) for various values of dimensionless unsteadiness S, the Deborah numbers β 1 and β 2 , when N t = 0.5, N b = 0.5, Pr = 7 and Le = 10.   It can be seen that when β 2 increases, both the temperature and concentration volume fraction distributions decrease thus diminishing the thicknesses of both the thermal and volume fraction boundary layers. Figs 12 and 13 show the effect of the thermophoresis parameter Nt on θ(η) and ϕ(η) for fixed S, β 1 , β 2 , Le and Nb. The temperature gradients in the boundary layer induces a thermophoretic force on the nanoparticles and that leads to a fast flow away from the stretching surface. Hence more fluid is heated away from the surface, and consequently, as N t increases, the temperature within the boundary layer increases. The fast flow from the stretching sheet carries with it nanoparticles leading to an increase in the mass volume fraction boundary layer thickness. It can also be observed that with an increase in the thermophoretic force, the nanoparticle fraction concentration profiles increase in the boundary layer before reducing to zero far from the surface. Here boundary layer separation occurs early at the stretching surface. Figs 14 and 15 show the influences of the random particle motion (represented by the parameter N b ) and the Lewis number Le on the nanoparticle volume fraction profiles. The nanoparticle Brownian motion at the molecular level plays a significant role in determining the thermal behaviour of the nanoparticle-fluid suspensions, Jang and Choi [45]. It is obvious that the nanoparticle volume fraction increases close to the stretching surface with increased Brownian motion and Lewis numbers before the boundary layer separation point. However, after the separation point the concentration volume fraction profiles decrease with an increase in both N b and Le.    Figs 16 and 17 show the effects of Deborah numbers β 1 and β 2 on skin friction and heat transfer coefficients. The skin friction coefficient increases with increasing β 1 leading to higher surface shear stresses, while increases in the skin-friction coefficient takes place with increasing β 2 . We observe that increasing β 1 leads to a decrease in −θ 0 (0) whereas increasing β 2 enhances the rate of heat transfer.
The variation of heat transfer coefficients with the thermophoresis parameter Nt is shown in Fig 18. It is clear that the thermal boundary layer thickness increases when the thermophoresis parameter Nt increases, and hence reducing the rate of heat transfer.
In order to apply the linear regression formula used by Kuznetsov and Nield [46], we used a set of 125 values of S, Nb and Nt with each S, Nb, Nt restricted to the space [0, 0.5] with a maximum error of less than 1%. Table 4 shows the linear regression coefficients and error bounds for the reduced Nusselt number. Here C s , C b , C t , are the coefficients in the linear regression estimate and ε is the maximum relative error defined by ε = j(Nu est −Nu)/Nuj. In this study, the minimum error occurs for small values of S, Nb and Nt.

Conclusions
In this paper we have studied the unsteady Oldroyd nanofluid flow over stretching surface. The classical boundary condition in which both the nanoparticle volume fraction and the   temperature are actively controlled has been substituted by the more realistic condition where the nanoparticle volume fraction is not controlled an the boundary. The effects of the governing parameters such as the unsteadiness parameter, the Deborah numbers in terms of relaxation and retardation times, the Prandtl number, the Brownian motion parameter, the thermophoresis parameter, the Lewis number on skin friction, heat transfer coefficients and fluid flow characteristics have been studied. Here β 1 represents the viscoelastic properties of the fluid and resists the motion of the fluid. The effects of the Brownian motion on the rate of heat transfer are negligible. The comparison between results obtained using the SRM and the QLM for skin friction and heat transfer coefficients showed a good agreement, with the SRM having converged at the sixth order up to six decimal places.  Table 4. Linear regression coefficients and error bounds for the reduced Nusselt number.