Stability analysis of triple solutions of Casson nanofluid past on a vertical exponentially stretching/shrinking sheet

The MHD two dimensional boundary layer flow of Casson nanofluid on an exponential stretching/shrinking sheet is considered with effects of radiation parameter, nanoparticles volume fractions (i.e. Fe3O4 and Ti6Al4V) and thermal convective boundary condition. The partial differential equations are transformed into ordinary differential equations by means of similarity transformations. The solutions of the transferred equations are achieved numerically with the help of shooting technique in Maple software. At different ranges of involved physical parameters, triple solutions are found. Therefore, stability analysis is performed by bvp4c in MATLAB to find the stable and physically reliable solution. Impacts of the physical parameter are presented through graphs and tables. Mainly, it is found that an increase in Casson and suction parameters decrease the corresponding velocity profiles while increase in Prandtl number, stretching/shrinking, and suction parameter decrease the temperature profiles. Furthermore, an increase in nanoparticles volumetric fraction, radiation and magnetic parameters as well as Biot number increase the temperature profiles and their thermal boundary layer thicknesses.


Introduction
The boundary layer flow and heat transfer characteristics of the fluids have wide range of applications in industrial and engineering zones. The processes involving flow on stretching and shrinking sheets (surfaces) have also achieved much importance due to its growing needs in industrial sectors for example, in glass production industries, paper production, polymer extrusion and plastic films, etc. In this regard, first significant attempt in study of stretching surface was reported by Crane. 1 The Crane's work was further extended by Chen 2 by adding mass transfer characteristics in various physical phenomena. The heat and the mass transfer characteristics on boundary layer fluid flow over exponentially stretching surface was examined by Magyari and Keller. 3 Elbashbeshy 4 considered the fluid flow past on the exponentially stretching surface. An investigation on two-dimensional steady boundary layer flow and heat transfer characteristics on exponential stretching surface with thermal radiation effect was considered by Bidin and Nazar. 5 Recently, few investigators have examined the problems of the fluid flows over stretching surface (see Khan et al. 6 and Zehra et al. 7 ).
Whereas, since many years, the heat transfer and thermal conductivity enhancement have become an important field of research. In this regard, most of the researchers have mainly focused on the study of nanofluids which are prepared by suspending nano-sized particles of different solids materials in different conventional fluids that are named as a base fluids of the nanofluids. These fluids can be utilized as a working fluids instead of conventional common fluids in future due to possessing high thermal conductivity, heat capacity, and heat transfer enhancement. In this regard, Choi and Eastman 8 were the first who prepared an innovative fluid by suspending solid nano-sized particles in the common base fluid. Lee et al. 9 stated that nanofluids possess outstanding heat transfer characteristics as compare to common fluids. Nanofluids can be used in various energetic systems of modern science and technology, like as natural convection enclosures, radiators, cooling of the nuclear systems, etc. Nanofluids are prepared by using Newtonian or non-Newtonian base fluids through suspensions of different types of solid nanoparticles that possess enhanced thermo physical properties as compare to base fluids. To understand the flow and heat transfer characteristics of nanofluids on mathematical viewpoint, Buongiorno 10 and Tiwari and Das 11 presented mathematical models. Furthermore, the boundary layer nanofluid flow on a nonlinear stretching sheet was studied by Rana and Bhargava. 12 Mabood et al. 13 examined the magnetohydrodynamic (MHD) nanofluid flow and examined the heat transfer characteristics for nonlinear stretching sheet. The studies of the boundary layer nanofluids flows over stretching/shrinking surfaces have been considered by many researchers (see Lund et al. 14 , Anwar et al. 15 , Nadeem et al., 16 and Abbas et al. 17 ).
In daily life, there are many applications related to fluids such as condensed milk, shampoos, tomatoes juice, printing ink, muds, and paints, etc. All these types of fluids show different types of the characters that are not understandable properly by the Newtonian theory of the fluids. Therefore, it is compulsory to study the non-Newtonians fluids to understand the phenomena of such type of fluids which are full with complexities. The properties related to non-Newtonian fluids may not be stated using a single model of non-Newtonian fluids properly. Therefore, different models have been introduced that are present in literature which are particularly categorized in three different models, namely integral, rate and differential type fluid. In 1959, a model proposed by Casson 18 that was later named as the Casson fluid flow model in which Casson fluid shows a yield stress. It shows an infinite viscosity at the zero shear stress while it is a type of the shear thinning fluid. In case of the Casson fluid, if applied shear stress is less than yield stress than no flow occurs, whereas, zero viscosity occurs at infinite shear rate, that is, when the shear stress is greater than the applied yield stress, fluid starts moving, whereas fluid show same behavior as a solid if the yield stress is greater than the applied shear stress. The purpose of the present study is to investigate the mixed convection magnetohydrodynamic (MHD) boundary layer flow and the heat transfer characteristics of Fe 3 O 4 and Ti 6 Al 4 V Casson based nanofluid on exponentially stretching/shrinking surface by using Tiwari and Das's 11 model. Present analysis is also modification of work done by Rohni et al. 19 and Hafidzuddin et al. 20 The effects of Casson fluid flow are discussed by many researchers (see Hafidzuddin et al. 20 , Amjad et al. 21 ).
In this article, similarity transformations are used to transfer partial differential equations into the system of ordinary differential equations that is requirement of the numerical method which is used here. At first, the obtained ordinary differential equations are numerically solved by employing shooting method in Maple software. The obtained results of the equations show that there exist three solutions for each value of the parameters which remains continue at different limits of the pertinent parameters. To check the reliability of the solutions, the stability analysis of the solutions is done by the bvp4c method in MATLAB. The obtained results of the present paper will prove supportive for researchers which are interested to find the impact of different physical fluid parameters on Casson-based nanofluids with the use of the different nanoparticles by applying Tiwari and Das's model. This study can be beneficial for those which are interested to work on nanofluids on mathematical and experimental basis because in present study two types of the nanoparticles are used and checked their flow and heat transfer characteristics. On other hand, in case of the rising multiple solutions, there has been given process of the stability analysis which indicates that among three obtained solutions, only one solution is reliable for the present study.

Formulation of problem
The two dimensional mixed convection magneto hydro-dynamic (MHD) boundary layer flow and heat transfer of Fe 3 O 4 and Ti 6 Al 4 V-Casson nanofluid over an exponential stretching/shrinking sheet with thermal radiations effect is considered. The Tiwari and Das's model is used to develop the equations of the present problem. A continuous magnetic field B(x) is used vertically across in the flow direction. The flow geometry and its coordinate system are denoted in the Figure 1. According to Amjad et al., 22 and Ali et al., 23 the Casson fluid flow rheological equations are considered as: where, t ij is the stress tension component, p is product of the deformation rate by itself, e ij is i, j ð Þ th component of deformation rate, p c denotes critical with value base on non-Newtonian model and the P y = m B ffiffiffiffi Under above assumptions, the continuity, momentum and the energy equations that defines the Casson nanofluid in accordance with Tiwari and Das model and following to Rohni et al. 19 can be written as: here, u and v represent the components of velocities along the x and y axes. The surface is assumed vertically stretched with velocity u w = U w e x L , here U w is constant of velocity and L is reference length. The surface temperature is taken as ity. Furthermore, l stands for stretching and shrinking parameter. While by Roseland approximation for the radiation, the radiative heat flux is written as: here, k Ã denotes coefficient of mean absorption while Stefan Boltzmann constant is denoted by symbol s Ã . The difference in temperature in flows are supposed to be sufficiently smaller so as T 4 should be stated as the linear function of T , applying the truncated form of Taylor series on ambient temperature T ' and then ignoring the higher order terms due to smaller values, we get: The associated boundary conditions are, The required similarity solutions are obtained with the help of similarity transformations which are In form of the velocity components, the stream function c is and magnetic field B(x) with constant B 0 for exponential surface is assumed to be With the help of the transformations given in equation 9 ð Þ, the system of the equations (4), 5 ð Þ, and (8) is written as the boundary conditions are Moreover, j = Gr Re x 2 is mixed convection parameter where Gr = is Biot number. The skin friction coefficient (C f ) and the local Nusselt number (Nu x ) both are considerable quantities which are take as, here, t w is the shear stress of wall and q w is heat flux from surface to fluid that can be expressed as and q w = k nf ∂T ∂y Using equations (9) and (14) in equations (13), we get, where, Re x = Lu w q f is a Reynolds number.

Stability analysis
There occurs triple solutions in this problem and in such cases, it is necessary to perform the stability analy- ∂T ∂t here, t is time. t is taken as a new dimensionless variable, by variable t, the equation (9) can be expressed as: where . By applying equation (18) in equations (16) and (17), it is obtained as: boundary conditions equation (12) can be written as, The perturbed basic solutions where, g shows the smallest eigenvalue, F h ð Þ and G h ð Þ both are small relative to f 0 h ð Þ and u 0 h ð Þ respectively. In the last, we get a system of the linearized eigenvalue problem, 1 Pr By setting t = 0, we get, 1 Pr with boundary conditions, The above linearized equations (25) and (26) with boundary conditions equation (27) are need to be solved. Furthermore, to find the smallest eigenvalues (g) by using the bvp4c solver in Matlab software. To obtain the smallest eigenvalues, named as g, one of the boundary condition is transferred in the from of the initial condition that was proposed by Harris et al. 24 Here, in this problem, It is mentioned that the smallest obtained negative values of g show the disturbance in initial growth, so the solution is said to be unstable. While, if obtained smallest value of g is positive, then the solution is said to be stable and physically feasible.

Numerical Method
The boundary value problem (BVPs) expressed in equations (10) and 11 ð Þ subjecting to initial and the boundary conditions specified in equation 12 ð Þ are solved by shooting technique which is explained by Meade et al. 25 A short procedure is defined as, along the conditions where, a 1 and a 2 are unknown initial conditions. Hence, the initial values of a 1 , a 2 is crucial that the solution might satisfy the boundary conditions F P h ð Þ ! 0, u h ð Þ ! 0ash ! ' of the present boundary value problem. The calculation is done by numerical technique named as shooting method in the Maple software. In present study, the boundary layer thickness (h ' ) is taken between 5 and 6 which provides a good convergence in solutions. The ranges of the velocity boundary layers thicknesses remained comparatively higher than thermal boundary layer thicknesses. Consequently, that provides three different values for f 0 0 0 ð Þ and Àu 0 0 ð Þ, respectively.

Results and discussion
The equations (10) and (11) that subject to the boundary conditions in equation (12) 20 which is presented in Table 1. The values of the thermo-physical properties of Fe 3 O 4 and Ti 6 Al 4 V nanoparticles as well as Casson fluid are presented in Table 2. Due to occurrence of the triple solutions, there has been done stability analysis to find the feasibility of the solutions that to be checked which one solution is stable and physically feasible. To obtain the stability analysis, we have used another technique that is three stage Lobatto III a formula. This formula is constructed in bvp4c by utilizing the finite difference code. Later on, the bvp4c solver functions are applied to find the results. The obtained smallest eigenvalues are presented in Table 3, which indicates the stability of first and physically feasible solutions due to having positive smallest eigenvalues. While, second and third solutions are found unstable that can not be physically feasible solutions due to having negative smallest eigenvalues.  Table 1. Comparative of the values of f 00 (0) and Àu 0 (0) with those of Hafidzuddin et al. 20 for different S when f = Rd = j = 0, l = 1, b ! ', Bi ! ', and Pr = 0:72.

S
Hafidzuddin et al. 19 Present       Figure 11 shows the variation in velocity profile of Fe 3 O 4 -Casson nanofluid by increase in the Casson parameter (b). In the first (stable) solution, an increase in b decreases velocity profile with its boundary layer thickness throughout the flow. While in second and the third solutions, velocity profiles are decreasing initially, but after a point they are increasing in both solutions. Actually, an increase in the b rises the dynamic viscosity. The yield stress that drags the nanofluid to the stretching surface which develops the resistance in flow due to rising the dynamic viscosity. In result, the velocity and the their boundary layer thicknesses are decreasing when parameter b is increased. Figures 7 to 12 show the effect of Casson parameter (b) on the temperature profile of      Figure 13 demonstrates that an increase in l, the velocity, and the boundary layers thickness are increased in first (stable) solution clearly. While in the second and in third solutions the fluctuations are observed where increasing and decreasing behaviors are observed. Figure 14 indicates the temperature profile and thermal boundary layer thickness decrease for larger values of l in first (stable) and third (unstable) solutions throughout the flow, while in second solution, it increases. The change in velocity and temperature profiles of Fe 3 O 4 -Casson nanofluid flow due to the suction parameter (S) is given in Figures 15 and 16, respectively. Figure 15 indicates the velocity profile and its boundary layer thickness decreases with increase in S in first (stable) solution. Actually, flow of the nanofluid comes nearer to solid sheet and velocity boundary layer becomes thinner when S is increased. While in second and third solutions, the velocity profiles decrease initially, but after a point they are increasing in both unstable solutions. It is noted from Figure 16, the temperature profile and the thermal boundary layers decrease when the rate of the suction is increased in all three obtained solutions. This is due to that the flow comes closer to boundary layer by increase in mass suction, thus,    thickness of velocity boundary layer takes the thinner form and also temperature losses easily with flow outward. Figure 17 indicates the relation of magnetic parameter (M) with the temperature profile of Fe 3 O 4 -Casson nanofluid flow. This figure demonstrates the velocity in first and third solutions increases throughout the flow with increase in magnetic parameter and same behaviours can be seen for the thickness of its momentum boundary layer in this particular case. Actually, an increase in magnetic parameter develops a resisting force (Lorentz force) in flow of the fluid, that acts in contrast to the flow direction. This resistance in flow becomes a main reason of enhancement in temperature profile when magnetic parameter is increased. While the second solution indicates the decreasing trend of the temperature profile with its boundary layer thickness. Figure 18 is drawn to express the effect of radiation parameter (Rd) on temperature profile. It can be seen from this figure, the temperature profile and thickness of its thermal boundary layers are increasing in first solution with increase in Rd. While in the second solution, temperature profile and its boundary layer thickness decrease throughout the flow and in third solution, it increases at the start but after a while it decreases as Rd is increased. The rise in radiation parameter emit heat energy toward the flow of the fluid, therefore, the temperature field is increased by heat energy. Thus, the high values of the radiation develops thickness of thermal boundary layer. The impact of Biot number (Bi) on temperature profile is illustrated in Figure 19. It is clear that the increasing value of the Bi increases the temperature profile in first solution and        Also, their velocity boundary layer thicknesses reduce with increase in both parameters. The higher values of the Prandtl number, stretching/shrinking, and suction parameter which decrease the temperature profile. An increase in nanoparticles volumetric fraction, radiation and magnetic parameters as well as Biot number rise the temperature profiles and related thermal boundary layer thicknesses.

Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding
The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: The authors extend their appreciation to the Deanship of Scientific Research at King Khalid University for funding this work through research groups program under Grant No. RGP.1/183/42.