BUOYANCY RATIO AND HEAT SOURCE EFFECTS ON MHD FLOW OVER AN INCLINED NON-LINEARLY STRETCHING SHEET

This paper numerically investigates the magnetohydrodynamic boundary layer convective flow of an electrically conducting fluid in the presence of buoyancy ratio, heat source, variable magnetic field and radiation over an inclined nonlinear stretching sheet under convective surface boundary conditions. The Rosseland approximation is adopted for thermal radiation effects and the non-uniform magnetic field applied in a transverse direction to the flow. The coupled nonlinear momentum, thermal and species concentration governing boundary layer equations are rendered into a system of third order momentum and second order energy and mass diffusion ordinary differential equations via similarity transformations with appropriate boundary conditions. The nonlinear, non-dimensional, well-posed boundary value problem is then solved with the implicit finite difference scheme known as Keller Box Method. A detailed study of the influence of the emerging dimensionless parameters governing the flow on velocity, temperature and concentration distributions is conducted. Also the evolution of skin friction coefficient, Nusselt number and Sherwood number values are depicted graphically. Numerical results are validated with some limiting cases documented in previously testified results, and good correlation is confirmed. This investigation is relevant to the delineation of space, astrophysical plasmas, polymer processing and extrusion of rubber and plastic sheets.


INTRODUCTION
The convective transport phenomenon in a porous medium received considerable attention stems from various engineering applications in geothermal reservoirs, thermal insulation engineering, petroleum industries, separation process in chemical industries and solar heating systems. Currently there is a significant activity in a special sub-category of boundary layer flows known as Sakiadis flows (1961) which deals with a wide variety of boundary layer problems with fluid flow over a stretching sheet in both Newtonian and non-Newtonian fluids. These flows are used in polymer synthesis, artificial fibers, hot rolling, paper production, chemical process, material processing, electrochemistry and many more possible areas. Crane (1970) studied the characteristics of boundary layer flow past a stretching plate. Flow arises due to stretching flat surface in three dimensional space have been studied by Wang (1984). Forced convection flows in stretching sheet with Marangoni effects was examined by Chen et al. (2007). Magyari and Chamkha (2008) studied MHD thermo solutal Marangoni convection in the presence of uniform magnetic field in the boundary layer approximation. Further, Chamkha (1997) adopted Finite Difference Method (FDM) to study Magnetohydrodynamic (MHD) free convective flow from a vertical plate with Hall effects embedded in a thermally stratified porous medium. It has been observed from earlier research that in fluid flow situations boundary circumstances are considered either at a particular specified surface heat flux or at a wall temperature. Hence, convective boundary condition effects arising in different boundary layer flow situations are important. Ramesh et al. (2015) discussed the effects of fluid particle suspension and convective boundary condition past a stretching sheet. Aljoufi and Ebaid (2016) elaborated the influence of convective boundary conditions on boundary layer slip flow over stretching sheet. For three dimensional non-Newtonian flow, Sulochana et al. (2016) examined similarity solution with convective boundary conditions. Ram Reddy and Pradeepa (2017) adopted spectral quasi linearization method (SQLM) to investigated nonlinear thermal convective flow of micropolar fluid under convective boundary conditions. Chakraborty et al. (2016) obtained numerical solution with convective boundary conditions for bio convection of a nanofluid flow. Extensive discussion of other applications recently studied by several researchers over stretching sheet which includes Sisko fluid flow employing Homotopy Analysis Method (HAM) Khan and Shahzad (2012), Second grade MHD fluid flow using Runge-Kutta (RK) sixth order integration scheme Das et al. (2016),  considered semi-infinite inclined and ideally transparent flat plate embedded in a porous medium with an Implicit finite difference scheme (IFDM), mixed convective couple stress fluid flow in a vertical channel with HAM by Kaladhar and Srinivasacharya (2014) and Prandtl-Eyring fluid flow with Keller-Box Method (KBM) Hussain et al. (2017). Further a significant amount of research on Newtonian and non-Newtonian fluid flows over stretching sheet were found in (Akbar et al., 2016;Kumar et al., 2017;Chamkha et al., 2010;Reddy and Gorla, 2017).
The most fundamental fluid mechanics problems widely found in nature are buoyancy driven flow due to release of heavy fluid into light fluid and also complexity arising from temperature dependent density and viscosity fields which in turn causes buoyancy. This buoyancy induced flow received considerable attention due to the potential applications such as heat exchanger design, petroleum production, chemical catalytic reactors and nuclear waste reactors and many possible areas. Soong (2001) analyzed the transport phenomena of developing thermal buoyancy effects in rotating non-isothermal flows. Utilizing nanofluids Khanafer et al. (2013) studied buoyancy driven heat transfer enhancement in two-dimensional enclosure. Rashidi et al. (2014b)

Frontiers in Heat and Mass Transfer
Available at www.ThermalFluidsCentral.org studied buoyancy effects on MHD flow of nanofluid (with Copper and Copper Oxide as nanoparticles) over a stretching sheet in the presence of thermal radiation by adopting fourth order RK based shooting technique. They concluded that Cu-water has more temperature than Cuo-water nanofluid while the buoyancy decrease the temperature of the nanofluid. Sheikholeslami et al. (2016a) studied the effect of MFD (Magnetic Field Dependent) viscosity on free convective nanofluid flow in an enclosure with bottom wall has constant heat flux element by using control volume based Finite Element Method (CVFEM) with linear triangular elements. Makinde et al. (2013) concluded that dual solution exists for shrinking case while studying stagnation point flow and heat transfer of a nanofluid past a convectively heated stretching/shrinking sheet with buoyancy effects. Anwar et al. (2012) studied the effects of buoyancy, solutal buoyancy and power law velocity parameters by adopting Keller-Box Method (KBM). Chamkha (2000) considered similarity equations governing the steady hydromagnetic boundary layer flow over an accelerating permeable surface ith buoyancy effects and then these equations are solved numerically with IFDM. Patrick and Paul (2010) for narrow vertical flat plate with uniform surface heat flux and plate edge conditions. Sheikholeslami et al. (2015) adopted CVFEM developed in FOTRAN with triangular elements to study the forced convection heat transfer in a lid driven semi annulus enclosure filled with Ferro nanofluid in the presence of non-uniform magnetic field and this extended to FHD (Ferrohydrodynamics) by the Sheikholeslami et al. (2016b), in this work they used CVFEM to study convection heat transfer in semi annulus under the influence of a variable magnetic field considering both FHD and MHD. They concluded that Kelvin force is more pronounced for high Reynolds numbers. Magyari and Chamkha (2013) reported an exact solution for the effects of buoyancy force and chemical reaction on micropolar fluid flow over a permeable stretched surface. Very recently Barletta et al. (2017) investigated unstable buoyant flow in a vertical porous layer considering convective boundary conditions and few related studies can be seen in Chemseddine et al. (2017) for molten PB-SN alloys.
Heat source/sink effects arise in many complex thermal technologies featuring high temperature differences. These include underground disposal of radioactive waste material, dissociating fluids in packed bed reactors, fire and combustion processes, blast furnaces, nuclear reactor fuel debris heat removal and in material fabrication of powders the working fluid heat generation and absorption effects are important. Early studies dealing with heat generation or absorption include the works of Chamkha (2003) Rashidi et al. (2011) used HAM to present the analytic solution to heat transfer of micropolar fluid through a porous medium with radiation effects and finally Chamkha et al. (2000) Hiemenz flow through porous medium. Rashidi et al. (2013) applied second law of thermodynamics to analyze an electrically conducting MHD nanofluid flow over a porous rotating disk in the presence of uniform magnetic field with RK fourth order based shooting technique and they proposed the path for optimizing the entropy. Thirupathi et al. (2017) adopted KBM to study effective Prandtl number and heat source effects on nonlinear stretching/shrinking sheet. The contribution made by several researchers for vertical surfaces are found in (Takhar et al., 2001;Kaladhar and Srinivasacharya, 2016;) however the study of inclined surfaces also important due to buoyancy forces that can be manipulated via inclination. Shit and Haldar (2011) investigated the action on thermal radiation Hall current over an inclined permeable stretching sheet. Recently Hasan et al. (2017) numerically investigated magnetohydrodynamic free convection flow past inclined stretching sheet using Nachtsheim-Swigert shooting technique with classical sixth order RK iteration. By using IFDM Takhar et al. (1999) studied the laminar boundary layer flow of an electrically conducting fluid past a semi-infinite flat plate with an aligned magnetic field effects as well as Takhar et al. (2002) solved MHD flow over a moving plate in a rotating fluid with IFDM. However, studies (Sandeep and Jagadeesh (2016), Rawi et al. (2016) and Amit (2017)) which focused on non-Newtonian flows along inclined surfaces demonstrate significant influence of inclination on MHD thermo fluid dynamic characteristics.
Although significant studies of magnetic, radiative and electrically conducting fluid flows from horizontal stretching sheets have been reported but in many metallurgical processes involve cooling of continuous strips or filaments and these strips are stretched when they are drawn through a quiescent fluid and also the final product quality improved by drawing such strips in an electrically conducting fluid subject to magnetic field. This motivates to consider the MHD boundary layer convective flow of an electrically conducting fluid in the presence of buoyancy ratio, heat source and radiation over an inclined nonlinear stretching sheet under convective surface boundary conditions. Therefore, the objective of present study is to extend the work with semi analytical method (HAM) reported by Rashidi et al. (2014a) and the numerical study (RKF-45 with scaling group transformation) reported by Ferdows et al. (2013). The system of coupled non-linear ordinary differential equations are solved by adopting an implicit finite difference scheme due to Cebeci and Bradshaw (1984). Extensive details of the mathematical formation of the problem, numerical method of solution, grid independence study and validation are presented in sections 2-4 respectively. The effects of various governing physical parameters on velocity, temperature and concentration profiles are discussed via figures and tables in section 5. Finally, a summary of noteworthy results is presented in section 6. The current simulations are relevant to paper production, glass fiber, polymer processing, cooling of metallic sheets in a metallic bath and many more.

MATHEMATICAL FORMULATION OF THE MODEL
Consider the steady, two-dimensional, incompressible, MHD, laminar and free convective flow of an electrically conducting fluid from an inclined permeable nonlinear stretching sheet, orientated at an angle α to the vertical. The physical model is illustrated in Fig.  1. The x-axis is directed along the continuous permeable stretching sheet and the y-axis is measured normal to the axis x − . The stretching of the sheet is induced by applying two equal and opposite forces simultaneously along the x-axis. By keeping the origin fixed, the sheet is stretched with nonlinear velocity , where c is a constant stands for characteristic stretching intensity.  1 It is assumed that flow takes place for 0 ≥ y and the permeable stretching sheet surface is heated by convention from a hot fluid with the temperature f T while the temperature of the electrically conducting ambient fluid takes the constant value ∞ T as ∞ → y and it is also assumed that the temperature varies in the axis x − direction also as is constant related to magnetic field. The external electrical field assumed to be zero and the electric field due to polarization of charges is negligible. The fluid properties are assumed to be constant except the density in the buoyancy terms in momentum equations and induced magnetic field strength is less when compared with applied external magnetic field hence this assumption is valid for the small Reynolds number. Finally, pressure gradient, Hall current, viscous and Ohmic dissipation are assumed to be neglected. Under these aforementioned assumptions with the Boussinesq approximations, the boundary layer governing equations following Rashidi et al. (2014) for Continuity, Momentum, Energy and species Concentration in the presence of thermal radiation, non-uniform magnetic field and heat generation can be shown to take the form: The following boundary conditions are imposed at the stretching sheet and in the freestream are: is the nonuniform magnetic parameter, g is the gravitation field, is thermal conductivity of the fluid and D is the coefficient of the mass diffusivity.
With an assumption for the optically thick boundary layer we adopt Rosseland's diffusion approximation for the radiative heat flux r q as given by y constant, * k is the Rosseland mean absorption coefficient. The term 4 T due to radiation can be expressed as a linear function of temperature itself. Therefore, 4 T can be approximated by Taylor's series about ∞ T after neglecting the higher order terms as Proceeding with the analysis, it is pertinent to introduce the following similarity transformations: Where the stream function ψ is defined as identically satisfies the continuity Eq. (1). Now by substituting Eq. (7) into the Eqs.
(2), (4), (6) and the corresponding boundary conditions of the problem considered in Eq. (5), the governing boundary layer equations reduce to the following system of highly coupled nonlinear ordinary differential equations, in terms of a single independent variable ) (η which are given by: Subject to the transformed boundary conditions for the stretching flow are given by: is the magnetic field parameter, Using Eq. (7) and Eq. (12) the skin friction coefficient, the local Nusselt number and local Sherwood number can be expressed as is the local Reynolds number.

NUMERICAL SOLUTION OF NONLINEAR BOUNDARY VALUE PROBLEM
The system of nonlinear coupled and inhomogeneous ordinary differential equations Eq. (8) -Eq. (10) subject to the boundary conditions in Eq. (11) are solved numerically using an implicit finite difference scheme known as the Keller-Box Method. Furthermore, this method is found to be suitable in dealing with nonlinear parabolic partial differential equations. Further details of this method are available for convection flows in the monograph of Cebeci and Bradshaw (1984). This method has four fundamental steps. The first step involves converting the Eqs. (8) -(10) into a system of first order ordinary differential equations. Thus, the coupled differential equations of third order in ) (η f and second order in ) (η θ and ) (η φ has been reduced to a system of seven simultaneous equations of first order for seven unknowns as follows. ; In the second step derivatives are approximated in the new system of first order equations with central difference approximations by considering the net rectangle in the plane x − η as shown in Fig. 2    and then cast into matrix vector form. Finally, the linearized algebraic equations are solved using a block tri-diagonal elimination scheme implemented in MATLAB software with the suitable initial solution. This method is unconditionally stable (Cebeci and Bradshaw (1984)), has a second order accuracy. For this iterative scheme to solve the system of equations, a convergence criterion is required. This is specified as follows: when the difference between two successive approximations is sufficiently small ( 7 10 − ≤ ) the solutions are taken to have converged to the requisite accuracy. Mathematical details are omitted for brevity.

GRID INDEPENDENCE STUDY AND VALIDATION OF RESULTS
In order to ensure the obtained numerical solutions are independent of mesh density i.e. grid specification, the effects of grid size on the are replaced by a sufficiently large value where the velocity, temperature and concentration profiles approach to zero. The MATLAB code is ran with three different grid sizes 801x801, 1001x1001 and 1201x1201 when 2 . 0 = S as shown in the Table 1 for stretching sheet and a very good agreement is observed between the profiles. Finally, in the present study, the far field boundary condition ∞ η is replaced by a sufficiently large value 10 max = η and for all the computations in this numerical study the adopted grid size is 1001x1001 with step size of 0.001. Furthermore, to ensure the correctness of the numerical results acquired through the implicit finite difference scheme, a comparison of skin friction coefficient and Nusselt number is made with the results obtained through semi analytical method (HAM) reported by Rashidi et al. (2014) Table 2. This allows verification of the Keller box solutions with those reported by Rashidi et al. (2014).
A comparison of skin friction coefficient, Nusselt number and Sherwood number is also made in the absence of heat source for electrically conducting MHD flow over vertical stretching sheet studied by Ferdows et al. (2013) wherein the scaling group transformation is employed for similarity and then RKF-45 with a shooting technique was employed to solve the resulting ordinary differential equations. These comparisons are presented quantitatively in Table 3. It is evident from Table 3 that present Keller box finite difference numerical values correlate closely with the solutions obtained by Ferdows et al. (2013) which confirms the accuracy of the present results and hence confidence in the present Keller-Box solutions is therefore justifiably high. Additionally in Table 4 Table 4 and 5 also provide a useful benchmark for other researchers who may wish to extend the present model and validate different numerical procedures against the present Keller box computations. . This circumstances is often encountered in engineering applications such as the cooling of electronic components and nuclear reactors and in many other areas.

RESULTS AND DISCUSSION
The parameter values to represent with physically realistic flows and in accordance with Rashidi et al. (2014) for which the numerical computations are carried out are presented in the respective figure legends. From all the figures it is also observed that for the case of suction ( 0 > S ) velocity, temperature and concentration profiles are more when compared with the profiles for injection ( 0 < S ) case. Furthermore, The CPU took 6.52 seconds to compute the velocity profiles, 5.37 seconds to compute the temperature profiles, 4.96 seconds to compute the concentration profiles for 1001 nodal points with the Intel core i3 processor under windows platform, which are computed by using the Matlab command tic; { Statements … … } toc;.

Velocity distributions
The

Temperature distributions
The variation of temperature distribution ( , high Biot number indicates higher internal thermal resistance of the sheet surface than the boundary layer resistance. As Bi increases, the temperature of the fluid enhances, although as per imposition of convective wall temperature boundary condition the greatest influence is at the sheet surface itself hence leads to increase in thermal boundary layer. The influence of variable magnetic field parameter M on hydromagnetic temperature distribution is represented in Fig. 8. It is inferred form Fig. 8 that the reduction in the momentum boundary layer thickness for larger values of M (=1.0, 1.5 and 2.0), the temperature distributions are enhanced. Further, it is worth mentioning here that the temperature distributions are starts with convective surface temperature and decays monotonically to asymptotic value and then satisfies the far field boundary condition. Fig. 9 shows the effects of several values of thermal radiation (Nr=0.5, 1.0 and 2.0) on temperature profiles. For optically thick fluids increase in Nr serves to energize the boundary layer and the input of thermal energy which is scaled with a cubic variation free stream temperature for thermal radiation compared with linear variation for thermal conduction. Therefore, an increase in Nr enhances the temperature magnitudes substantially due to increase in the thermal radiative heat transfer.  Fig. 10. From figure it is evident that an increase in Pr , induces a significant reduction in the temperature and therefore, cools the magnetic fluid regime, thereby decreasing thermal boundary layer thickness. Further, greater Pr corresponds to lower thermal conductivity of the fluid therefore as Pr increases leads to suppression in the thermal conduction heat transfer and hence plummet in temperature. The graphical representation of dimensionless ) (η θ for some representative values of steady state surface dependent heat source (Q) parameter is shown in Fig. 11. An increase in Q (=0.5, 1.0 and 1.5) enhance the temperature due to the fact that the addition of Q introduces thermal energy in to the flow, which then energizes the boundary layer and elevates the thickness of thermal boundary layer. The effect of Sc on temperature is shown in the Fig. 12. With increase in Sc also enhances the temperature. Physically Sc behavior on temperature has opposite with the behavior of Pr on temperature since the Schmidt number embodies the ratio of momentum to mass diffusivity. Therefore whenever Sc is increased, the thickness of thermal boundary layer also increases.  Fig. 13. Whenever M increased the concentration profiles are also increased due to reduction in the thickness of the momentum boundary layer.

Concentration distributions
The Schmidt number ( Sc ) represents the relative effectiveness of momentum diffusion in velocity boundary layer to species diffusion in concentration boundary layer (

CONCLUSIONS
Magnetohydrodynamic boundary layer convective flow of an electrically conducting fluid in the presence of buoyancy ratio, heat source, and variable magnetic field and radiation over an inclined stretching sheet under convective surface boundary conditions has been studied theoretically and numerically in this paper. The governing coupled nonlinear momentum, thermal boundary layer and species concentration equations are transformed into coupled ordinary differential equations by using similarity transformations. The Keller-Box implicit finite difference method is adopted to solve the emerging nonlinear boundary value problem. The originality of the present study entails an elaboration of the influence of governing parameters on velocity, temperature and concentration distributions for suction (porous case) and injection (nonporous) cases with variable magnetic field for inclined stretching sheet under convective boundary conditions. Verification of the Keller box solutions has been achieved via comparison with previously published reports and a very close correlation is observed. The principal findings of the current simulation may be summarized thus: • The velocity profiles are decreased with increasing values of magnetic field parameter, inclination angle parameter and Schmidt number and also for increasing values of buoyancy parameter velocity distributions are also increased for both suction and injection cases. • The temperature distribution values are decreased with Prandtl number where as an increasing in Biot number, Magnetic field parameter, Schmidt number, radiation parameter and heat source the temperature values are increased for both suction and injection cases. • The concentration distribution values are decreased with Magnetic field parameter where as an increasing in Schmidt number the concentration values are increased for both suction and injection cases. It is also observed that for the case of suction velocity, temperature and concentration distributions are more pronounced than the injection case. • The skin friction coefficient and Nusselt number are increased by increasing the Biot number and heat source parameter. However Sherwood number is decreased with increasing values of the Biot number and heat source parameter.