Unsteady Formulations for Stagnation Point Flow Towards a Stretching and Shrinking Sheet with Prescribed Surface Heat Flux and Viscous Dissipation

The unsteady stagnation point flow and heat transfer with prescribed flux towards a stretching and shrinking sheet with viscous dissipation is studied. Similarity transformation is adopted to initially convert the governing differential equations into nonlinear ordinary differential equations. The two-point boundary value ordinary differential equations (ODE) are subsequently converted into partial differential equations by introducing a time-marching scheme. A Crank-Nicolson NewtonRichtmeyer scheme is employed to discretize the resulting equations. Initial guesses are made for the dependent variables and the solution advanced in time until temporal variations of the scalar profile are diminished and the steady-state solutions satisfy the similarity equations. A variation of the heat flux at one of the boundaries produced noticeable variations in the temperature field that can be related to the magnitude of the Prandtl number and velocity ratio parameter.


Introduction
The impact of a viscous fluid on a solid object (Fig. 1) constitutes a prototypical stagnation point flow. Stagnation point flow and heat transfer with prescribed flux have continued to receive considerable attention over the years because of their wide application. The heat transfer over a surface becomes quite interesting when the heat flux boundary condition is taken into account particularly for those problems involving the cooling of electrical and nuclear components, space shuttle re-entry into the earth's atmosphere, melt-spinning processes etc. This class of problems often encounters rapid heat-flux changes, and possible meltdown hence overheating and burnout constitute very important design considerations. As a result, accurate prediction of the temperature profile becomes vital (Bejan [1]).
The study of viscous flow near a solid object can be traced back to Hiemenz [2] similarity analysis. That is why the stagnation point flow is sometimes referred to as Hiemenz flow. He reduced the Navier-Stokes equations governing the flow to third order ordinary differential equations by using similarity transformations. An extension of his work was carried out by Homann [3]. Hiemenz flow was further studied by Na [4] using the method of finite differences. Due to its fundamental nature and in combination with its practical importance, the classical problem has been generalized in several ways to include diverse physical effects. These comprise viscous or inviscid, forward or reverse flows over surfaces. For example, the study of flow over a stretching sheet has several industrial applications especially in paper and glass fiber production, and in the extraction of polymer sheets drawn through a stagnant fluid with a controlled cooling system. One of the first attempts to look into the flow of incompressible fluid over a linearly stretching sheet was conducted by Crane [5] who provided similarity solutions for the case of a steady two-dimensional incompressible boundary layer flow produced by a stretching sheet that moves in its own plane. The heat transfer component of this study was carried out by Carragher and Crane [6], while the mass transfer effect including suction and blowing was investigated by Gupta and Gupta [7]. Chiam [8] studied the stagnation point flow towards a stretching sheet when the stretching velocity is identical to the free stream velocity. In a similar work, Sakiadis [9,10] considered the case where the stretched surface was assumed to move with uniform velocity while Mahapatra and Gupta [11] found out that the profile of the boundary layer is a function of the stretching sheet parameters and the angle of incidence. Several other authors can be credited for work in this area. Dulal and Hiremath [12] studied the heat transfer characteristics in a laminar boundary layer flow of an incompressible viscous fluid over an unsteady stretching sheet placed in a porous medium in the presence of viscous dissipation and internal absorption or generation of heat. Relevant references in this area can be found in the books by Ingham and Pop [13], Nield and Bejan [14]. Relating heat transfer to fluid flow has provided a link between the application of fluid numerical simulation to the industry, This can be observed in such diverse areas as navigation, aviation, rocket propulsion, space exploration, ship building, food processing, and auto manufacture to mention just a few. The heat transfer characterisitcs of a stretching sheet becomes significantly important when the heat flux boundary condition is taken into consideration. Dutta et al. [15] examined the effects of uniform heat flux on the temperature field due to stretching. A similar study was carried out by Alli and Magyari [16] for power-law velocity and temperature distributions. The effect of a magnetohydrodynamic (MHD) flow on a stretching sheet was studied by Ishak et al. [17]. Similar work was carried out by Anderson [18] his work was later extended with the inclusion of mass transfer by Char [19]. Mahapatra and Nandy [20] provided a momentum and heat transfer solution involving a shrinking sheet for the case of an MHD axisymmetric stagnation point flow. Further work on hydromagnetic flows includes those of Makinde and Charles [21], Makonde [22]. Conjugate convective flows apply in cases involving internal heat generation. These often lead to obtaining useful design parameters [23, 24, and 25].
The aim of this study is to numerically investigate the problem of an unsteady stagnation point flow towards a stretching and shrinking sheet with prescribed surface heat flux and viscous dissipation. The governing differential equation is obtained from the vorticity transport equation and converted into its third order stream function analog. By using a combination of analytic techniques, nondimensionalization and appropriate boundary conditions it was possible to arrive at a second order nonlinear transient partial differential equation. The resulting system was solved to steady state by an implicit finite difference Crank-Nicolson-Newton-Richtmeyer scheme. The effects of various non-dimensionless parameters on the flow and heat transfer profiles are presented graphically and discussed.

Mathematical Formulations
We seek the governing equation for boundary layer flow towards a stagnation point flow (Fig. 1). The x and y axes are taken along the normal and along a shrinking stretching surface. It is assumed that the free stream velocity is given by x y y Boundary Conditions: Two boundary conditions on follow from the no-slip boundary constraint at the surface. The other two concur with the approach of the viscous-flow solution to inviscid flow far enough from the surface.
ψ is a constant in the inviscid flow solution. Because the first derivative of the streamfunction approaches a constant value 'A' asymptotically, we can deduce that its second derivative approaches zero as y ( ) , x y ψ ψ approaches infinity. This observation follows the logic that once a derivative of a function approaches a constant value, all higher order derivatives of that function will approach zero asymptotically. We make a further assumption that the streamfunction comprises both the viscous and inviscid components of flow and as a result; it can be expressed as Substituting this assumed form of ψ into equation (1) With accompanying boundary conditions: We adopt Wilcox [26] approach, and introduce the following dimensionless independent and dependent variables ( ) ( ) (3) and (4) to produce: with the following boundary conditions Equation (5) is the so called Hiemenz equation. A slightly different approach has been adopted for its derivation. One of the reasons for this is the relative ease with which boundary conditions are handled as well as the avoidance of an iterative procedure to arrive at an appropriate initial condition. Equation (5) facilitates the computation of the velocity field, but for many computations, we also need the pressure field. To initiate this, we make another assumption that the pressure comprises the static and dynamic components and can be described as follows: For the classical case of stagnation pressure, at 0 x y = = , ( ) 0 0 y F y = = , and the dynamic component is completely obviated. But we relax the conditions a little bit by considering the fact that as with velocity we expect the pressure to approach the inviscid flow solution in the limit as y → ∞ , given this criterion, the pressure should satisfy the condition: Since a clear concept of the pressure dynamics has been determined, we need an explicit equation for pressure to validate our assumptions. Taking the divergence of the momentum in a two-dimensional form of the Navier-Stokes equation and after rearranging terms, we obtain the Poisson equation for pressure Though it can be assumed that once the velocity field has been calculated from the streamfunction values, equation (9) can be solved to yield the pressure, nevertheless, it still has to be put in a form that satisfies the physics of the flow as dictated by equations (7) and (8). Equation (7) is further expanded as follows: and the function ( ) Y y satisfies: with the following boundary conditions: Equation (11) guarantees that 0 p be the true stagnation pressure while the pressure gradient matches that of the inviscid flow in a region far above the surface. In order to complete our solution of the pressure equation, we need to recast equation (10b) into its dimensionless form. We define: and substituting into equation (10b) we obtain: accompanied by the following boundary conditions Using equation (5), the squared gradient of the dimensionless streamfunction in equation (12) is eliminated. Integrating twice and using equation (13), we obtain a much simplified form of the pressure profile: We note that the original equation (equation (1)) for the streamfunction was fourth order and represents the viscous analog of the Laplace's equation for the streamfunction. In addition it satisfies the continuity equation automatically. After integrating once, equation (1) was reduced to third order ODE (equation 5). The overriding approach found in literature is to reduce the governing third order ODE to systems of lower order Odes before applying numerical techniques like the Runge Kutta method. We adopt a different approach here by reducing equation (5) to its second order analog before adopting a time marching procedure. We recall that the dimensionless streamfunction is related to the ( ) ' .
horizontal velocity u f η = This facilitates the formulation of a time marching scheme with the boundary conditions: For an initial condition, a guess is made by assuming a quadratic profile for the velocity profile near the surface of the plate Next we consider the boundary layer equation for heat transfer where ρ is the density, p c is the specific heat, µ is the kinematic viscosity, T is the temperature, x, y are the independent variables. The boundary conditions for the velocity and temperature fields are To solve the energy equation we obtain the velocity field from the solution of the momentum equation. These are then substituted into equation (16) in terms of the similarity variables and f η . We initiate this process by converting equation (16) and (17) into dimensionless forms by introducing: to finally obtain a dimensionless energy equation together with the boundary conditions:

Discretization
For this study, we have elected to convert equation (18) to a time marching numerical scheme by introducing a temporal derivative.
We employ the Crank-Nicolson Newton-Richtmeyer scheme for discretization, whose general form for a onedimensional partial differential equation for nonlinear diffusion, convection and source term is given as: Using a forward difference for the time derivative: Without any loss in generality; assume ( ) ϑ φ φ = , then equation (20) can be written as: The diffusion term can now be obtained as follows: By the same token Following this approach, a nonlinear reaction term is discretized as: where ϒ is the reaction coefficient, and χ is the power to which the nonlinear term is raised. We note that is now the new dependent variable and is related to the desired dependent variable by 1 1 constitute the tridiagonal coefficient matrix and together with the known right hand side known vector are solved with the Thomas algorithm for the dependent variables until steady state is achieved. The discretization of the momentum and energy equations follow the stencil set out in equation (20) and are much easier to apply for linear terms.

Results and Discussions
We  the pressure profile satisfies the boundary condition as specified in equation (11). The flow results are in agreement with those in page 527 Wilcox [26]. Since this work is heat flux based, we shall be focusing a lot of on temperature and temperature gradient profiles for different parameter values. The effects of Prandtl number, negative and positive velocity ratios for an isothermal left side boundary condition on the temperature and temperature gradient profiles are illustrated in Figs (3-6). For Fig. 3, each temperature profile, initially decreases and then ultimately goes to zero for large values of η . The temperature gradient profiles for Fig. 4 exhibit an initial decrease in gradient, then a rise followed by an asymptotic approach to zero for all the profiles. Based on these observations, it is worthwhile to note that considering the values of Prandtl numbers, if the momentum diffusivity is high (i.e. free movement of fluid due to density differences) and the thermal diffusivity is low (i.e. for low thermal conductivity fluids), the overall heat transfer is facilitated by natural convection. Prandtl number is relatively high for such cases. This is what obtains for natural convection of certain fluids, for example water.    The impact of a positive velocity ratio on the temperature and temperature gradient profiles can be seen in Fig. 5 and Fig. 6. As expected, both the temperature and the temperature gradient decrease along the surface in accordance with the specified boundary condition. But the positive velocity ratio effect on the zero flux boundary condition is profound. This is as a result of a significant drop in temperature and the alteration of the shape of the boundary layer. Hence an increase in the velocity ratio causes the temperature of the fluid to decrease. For such a case, the thickness of the thermal boundary layer decreases due to the fluid at the surface having a temperature lower than that of the flat plate.   also more prominent. With regards to the magnitudes of the scalar profiles and the assigned values of the Prandtl number, we note that if the momentum diffusivity is low (i.e. there is not enough density differences to enhance the movement of the fluid) and the thermal diffusivity is low (i.e. for fluids with low values of thermal diffusivity), then the transfer of heat from the from the boundary will not be facilitated as the fluid is not only incapable of absorbing heat but also cannot move freely due to insufficient momentum from convection currents. The profiles we encounter here are better assessed if we compare Figs. 7 and 3 to see how the introduction of a positive flux at the left boundary (with all the other parameters remaining the same) results in an overall relative increase in temperature within this region. It is obvious that we encounter a vector quantity here whose direction informs us that heat is drawn out of the plate.
Figs 9 and 10 show the effect of a positive velocity ratio on the temperature and flux profiles. An increase in velocity ratio results in a significant temperature loss implying that the thermal boundary layer decreases. Again it should be remembered that because of its dependence on the fluid velocity, the temperature field is affected by the velocity ratio. As the velocity ratio increases, the temperature of the fluid decreases. This is because the positive value of the velocity ratio promotes heat transfer to the fluid. Fig. 10 shows a monotonic rise in temperature gradient until it gets to a point where the temperatures of all the fluids become zero in agreement with the specified Dirichlet boundary condition at the right end of the plate.   Fig. 7 records a much higher temperature than Fig. 11. However Fig. 11 displays much steeper slopes than Fig. 7. Fig. 12 gives values of the flux and reflects the heat transfer trend of Fig. 11. Fig. 13 displays the temperature profiles for a unit negative flux specification at the left boundary, keeping the magnitude of the velocity ratio the same. This gives us a handle on the magnitude of heat flow as well as the direction. A comparison made between Figs 13 and 9 yields valuable information concerning the link between the temperature profiles displayed in each case as well as magnitude and signs (directions) of the flux specification. Figs 10 and 14 give further information concerning the flux. While Fig. 9 shows a monotonic decrease in temperature for all values of Prandtl number, Fig. 13 displays an exactly opposite trend. Each of these cases is in consonance with the vector interpretation of the Fourier's law. However, one thing they both have in common is that at the region 3 η > all the profiles move towards satisfying the Dirichlet boundary condition specification at the right boundary.

Conclusions
In this work, we have presented unsteady numerical formulations for stagnation point flow towards a sheet with prescribed surface heat flux and viscous dissipation. Since more emphasis is devoted to the energy component of the formulation, several graphical plots are made to illustrate the temperature and the temperature gradient fields. The numerical results obtained are graphically displayed and discussed for various problem parameters considered in the problem formulation. Numerical experiments involving flux variations are made on the left of the problem domain. For each of these results, we can observe that as we move far from the sheet, the effects of viscous dissipation vanish and the profiles intersect each other at some point.
Some of the basic features of this model can be applied to study a variety of scenarios involving more complex formulations such as heat characteristics for Newtonian and non-Newtonian fluids, flows in porous media and MHD slip flow to mention just a few.