Triple Solutions and Stability Analysis of Micropolar Fluid Flow on an Exponentially Shrinking Surface

: In this article, we reconsidered the problem of Aurangzaib et al., and reproduced the results for triple solutions. The system of governing equations has been transformed into the system of non-linear ordinary di ﬀ erential equations (ODEs) by using exponential similarity transformation. The system of ODEs is reduced to initial value problems (IVPs) by employing the shooting method before solving IVPs by the Runge Kutta method. The results reveal that there are ranges of multiple solutions, triple solutions, and a single solution. However, Aurangzaib et al., only found dual solutions. The e ﬀ ect of the micropolar parameter, suction parameter, and Prandtl number on velocity, angular velocity, and temperature proﬁles have been taken into account. Stability analysis of triple solutions is performed and found that a physically possible stable solution is the ﬁrst one, while all leftover solutions are not stable and cannot be experimentally seen.


Introduction
Generally, the investigation of the non-Newtonian fluid flow in two-dimensional problems is a hard task when multiple solutions are attempted to find because equations are related to high nonlinear terms. In spite of these challenges, the researchers are making efforts to tackle these problems for multiple solutions due to their wide range of applications in various science and industrial fields. The category of non-Newtonian fluids which pacts with the suspended micro-rotational particles, is known as the micropolar fluid. Eringen introduced the theory of micropolar [1,2]. He explained the impacts of couple stresses and local rotational inertia that cannot be described by the standard equations of Navier-Stokes. The micropolar equations are mathematically described for the theory of porous media and the theory of lubrication in the books by Lukaszewicz [3] and Eringen [4]. There are various applications of the micropolar fluids, for example, liquid crystals, particle suspension, animal blood, lubrication, turbulent shear flows, and paints. Lok et al., [5] considered the stagnation point flow of micropolar nanofluid and then succeeded to find the dual solution. It is also stated that only the noticed. This paper is divided into six sections; Section 1 is for the brief introduction of micropolar fluid and multiple solutions. Mathematical formulation, derivation of stability analysis, and methodology are kept in Section 2, Section 3, and Section 4 respectively. Section 5 is for the result and discussion and Section 6 is for the conclusion.

Mathematical Formulation
An incompressible laminar boundary layer two-dimensional flow of the micropolar fluid over the exponentially shrinking sheet has been considered. The corresponding velocities of x and y-axes are u and v. The shrinking velocity is assumed to be u w (x) = −U w e 2x . The temperature of the sheet is taken to be T w (x) = T ∞ + T 0 e x 2 , as shown in Figure 1. The N = N(x,y) is supposed as the angular velocity. The respective boundary layer movement equation, along with micro rotations and the heat transfer equations can be expressed as vectors in accordance with the abovementioned assumptions.
, 0], the micro-rotation vector is N, ρ stands for fluid density, µ for viscosity coefficient, κ for vertex viscosity, j is the density of micro-rotation, and γ stands for micropolar constant. We get following boundary layer equations according to the scale analysis.
Subject to these boundary conditions Now, we look for similarity transformation variables in order to transform Equations (6)-(8) with boundary conditions (9) Crystals 2020, 10, 283 4 of 14 By applying Equation (10) in Equations (6)-(9), we have the following system of similarity transformed ordinary differential equations 1 Pr Subject to boundary conditions where prime denotes the differentiation with respect to η, the micropolar material parameter is K = κ µ , Prandtl number is Pr = ϑ α , and suction is . The physical quantities of interest include skin friction, the stress of local couples, and the local number of Nusselt, which are described as By applying similarity transformation variables (10) in Equation (15), we have where Re x = xu w /ϑ is the local Reynolds number.
Crystals 2020, 10, x FOR PEER REVIEW 4 of 14 By applying Equation (10) in Equations (6)-(9), we have the following system of similarity transformed ordinary differential equations Subject to boundary conditions where prime denotes the differentiation with respect to , the micropolar material parameter is = , Prandtl number is = , and suction is = − ⁄ . The physical quantities of interest include skin friction, the stress of local couples, and the local number of Nusselt, which are described as By applying similarity transformation variables (10) in Equation (15), we have where = ⁄ is the local Reynolds number.

Stability Analysis
According to Nasir et al., [28] and Rana et al., [25], we need to introduce the unsteady form of Equations (6)-(8) in order to perform stability test,

Numerical Methods
The governing ODEs are highly non-linear and, therefore, we adopt the numerical approach in order to solve Equations (11)- (14) and Equations (26)- (19). In this study, two methods have been employed, namely shooting method and Three-stage Lobatto III-A formula, which were used in many research articles of same authors previously [refer to 17,20,23,26]. The descriptions regarding these methods are explained below.

Shooting Method
The shooting technique along with the Runge Kutta method of the fourth order is employed in order to obtain the numerical solutions of Equations (11)-(13) subject to the boundary conditions. Shooting method helps to reduce the third order ODEs (11)-(13) into the first-order ODEs, such that with conditions where α 1 , α 2 , and α 3 are called as unknown initial conditions. These three missing values α 1 , α 2 , and α 3 have to be obtained by using different shoots; this process of shoots will be continue until the profiles of the f (η) → 0; h(η) → 0; and θ(η) → 0 are satisfied the boundary condition η → ∞ . Maple (18) software has been used to convert the system of the third order ODEs into the system of the first order ODEs, for this process shootlib function is built-in Maple. Using RK method solves the system of the first order ODEs. Further, a detailed discussion about the shooting method with Maple software can be seen in the paper of Meade et al. [34].

Three-Stage Lobatto III-A Formula
Three-stage Lobatto III-A formula is built in BVP4C function with aid of C 1 piece-wise cubic polynomial in the finite difference code. According to Lund et al., [35] and Raza et al., [36], "this collocation polynomial and formula offers a C 1 continuous solution in which mesh error control and selection are created on the residual of the continuous solution. The tolerance of relative error is fixed 10 −5 for the current problem. The suitable mesh determination is created and returned in the field sol.x. The bvp4c returns solution, called as sol.y., as a construction. In any case, values of the solution are gotten from the array named sol.y relating to the field sol.x". In addition, Figure 2 explains the algorithm of the method for stability analysis of the solutions.

Result and Discussion
Aurangzaib et al., [29] solved these Equations (11)- (14) and found dual solutions, which was the main contribution of authors. Aurangzaib et al., [29] have given some strong statements during the investigation of critical points. Before going to detail, we check the accuracy of our method (shooting method) by comparing our results with previously published literature in Table 1 and found excellent agreement with them. Moreover, we also qualitatively compared our results with Aurangzaib et al., [29] in Figure 3 and found in the good agreement. This gives trust in our numerical calculation and urges us to additionally contemplate this problem. Aurangzaib et al., [29] found dual solutions and stated in the result and discussion section that "the present study shows that for K = 0.1, i.e., for micropolar fluid, the similarity solutions exist when S ≥ 2.3231 and no similarity solution exists for S < 2.3231". Firstly, we would like to clarify that there exist triple solutions not dual solutions; secondly, there is a range of multiple solutions and single solution. From Figure 4, a conclusion can be made that there exists multiple similarity solutions when S ≥ 2.3224 and only a single similarity solution exists for S < 2.3224 when K = 0.1. However, the range of multiple similarity solution is 2.3769 ≤ S and there also exists only a single solution when S < 2.3769 when K = 0.2 (see Figure 4). Furthermore, skin friction increases as suction is increased in the third solution. It is worth mentioning here that there is no range of no solution. This is one of the big reasons that insist us to reconsider and re-examine the whole problem by reproducing all of the results because there are triple solutions in order to provide true knowledge to the readers and researchers.

Result and Discussion
Aurangzaib et al., [29] solved these Equations (11)- (14) and found dual solutions, which was the main contribution of authors. Aurangzaib et al., [29] have given some strong statements during the investigation of critical points. Before going to detail, we check the accuracy of our method (shooting method) by comparing our results with previously published literature in Table 1 and found excellent agreement with them. Moreover, we also qualitatively compared our results with Aurangzaib et al., [29] in Figure 3 and found in the good agreement. This gives trust in our numerical calculation and urges us to additionally contemplate this problem. Aurangzaib et al., [29] found dual solutions and stated in the result and discussion section that "the present study shows that for K = 0.1, i.e., for micropolar fluid, the similarity solutions exist when S ≥ 2.3231 and no similarity solution exists for S < 2.3231". Firstly, we would like to clarify that there exist triple solutions not dual solutions; secondly, there is a range of multiple solutions and single solution. From Figure 4, a conclusion can be made that there exists multiple similarity solutions when S ≥ 2.3224 and only a single similarity solution exists for S < 2.3224 when K = 0.1. However, the range of multiple similarity solution is 2.3769 ≤ S and there also exists only a single solution when S < 2.3769 when K = 0.2 (see Figure 4). Furthermore, skin friction increases as suction is increased in the third solution. It is worth mentioning here that there is no range of no solution. This is one of the big reasons that insist us to reconsider and re-examine the whole problem by reproducing all of the results because there are triple solutions in order to provide true knowledge to the readers and researchers.        Figure 5 illustrates the effect of micropolar parameter on the h (0). The effect of local couple stress enhanced as the suction increases in the first and third solutions because increasing suction creates additional resistance in the flowing fluid inside the boundary layer. However, increments in the material parameter produce more coupling of stress. Figure 6 shows the nature of the heat transfer rate for various values of the suction. It has been examined that the heat transfer rate increases in the first and second solutions for the higher values of the suction parameter, while the third solution shows opposite compliance.  Figure 5 illustrates the effect of micropolar parameter on the ℎ′(0). The effect of local couple stress enhanced as the suction increases in the first and third solutions because increasing suction creates additional resistance in the flowing fluid inside the boundary layer. However, increments in the material parameter produce more coupling of stress. Figure 6 shows the nature of the heat transfer rate for various values of the suction. It has been examined that the heat transfer rate increases in the first and second solutions for the higher values of the suction parameter, while the third solution shows opposite compliance.     Figure 5 illustrates the effect of micropolar parameter on the ℎ′(0). The effect of local couple stress enhanced as the suction increases in the first and third solutions because increasing suction creates additional resistance in the flowing fluid inside the boundary layer. However, increments in the material parameter produce more coupling of stress. Figure 6 shows the nature of the heat transfer rate for various values of the suction. It has been examined that the heat transfer rate increases in the first and second solutions for the higher values of the suction parameter, while the third solution shows opposite compliance.   Finally, we plot Figures 7-10 to show the existence of triple solutions of velocity, microrotation, and temperature profiles for different values of material parameter K. In Figure 7, the dual nature of behavior has been noticed in the first solution. Velocity profile increases in the second solution when K is increased; the physical material parameter reduces the effect of drag force due to that thickness of the momentum boundary layer enhanced. On the other hand, the opposite trend has been observed in the third solution. Dual behavior can be noticed microrotation profile in Figure 8 for all solutions. The thickness of thermal boundary layer increases in the first and second solutions as material parameter K is increased, as in Figure 9, since the non-Newtonian parameter produces more viscosity, decreases the velocity of profiles, and forces fluid flow to stay on the hotter surface, as a result temperature of fluid increases and the boundary layer becomes thicker. However, the opposite behavior is noticed for the third solution. The Prandtl effect on the temperature distribution is depicted in Figure 10. It is observed that the temperature of fluid diminishes for the higher values of the Prandtl number for all solutions. Physically, it can be explained, as the Prandtl number (Pr = µc p k ) has an inverse relationship with the thermal conductivity and, consequently, diminishes the thickness of the thermal boundary layer.
Finally, we plot Figures 7-10 to show the existence of triple solutions of velocity, microrotation, and temperature profiles for different values of material parameter K. In Figure 7, the dual nature of behavior has been noticed in the first solution. Velocity profile increases in the second solution when K is increased; the physical material parameter reduces the effect of drag force due to that thickness of the momentum boundary layer enhanced. On the other hand, the opposite trend has been observed in the third solution. Dual behavior can be noticed microrotation profile in Figure 8 for all solutions. The thickness of thermal boundary layer increases in the first and second solutions as material parameter K is increased, as in Figure 9, since the non-Newtonian parameter produces more viscosity, decreases the velocity of profiles, and forces fluid flow to stay on the hotter surface, as a result temperature of fluid increases and the boundary layer becomes thicker. However, the opposite behavior is noticed for the third solution. The Prandtl effect on the temperature distribution is depicted in Figure 10. It is observed that the temperature of fluid diminishes for the higher values of the Prandtl number for all solutions. Physically, it can be explained, as the Prandtl number ( = ) has an inverse relationship with the thermal conductivity and, consequently, diminishes the thickness of the thermal boundary layer.       Table 2 shows the smallest eigenvalues for the selected values of S and K. The positive smallest eigenvalue makes the initial disturbance decay and, in this way, the flow becomes stable. Conversely, the negative smallest eigenvalue outcomes in an initial growth of disturbance, in this manner, the flow is unstable. It is seen from Table 2 that is negative for the second and the third solutions, while positive for the first solution. Thus, the second and the third solution are not stable, and the first solution is stable. From this discussion, it can be concluded that the first solution of Aurangzaib et al., [29] is not stable and not physically realizable; therefore, in this stage, it could be said that the first solution is actually the second or the third solution.   Table 2 shows the smallest eigenvalues ε for the selected values of S and K. The positive smallest eigenvalue makes the initial disturbance decay and, in this way, the flow becomes stable. Conversely, the negative smallest eigenvalue outcomes in an initial growth of disturbance, in this manner, the flow is unstable. It is seen from Table 2 that ε is negative for the second and the third solutions, while positive for the first solution. Thus, the second and the third solution are not stable, and the first solution is stable. From this discussion, it can be concluded that the first solution of Aurangzaib et al., [29] is not stable and not physically realizable; therefore, in this stage, it could be said that the first solution is actually the second or the third solution.

Conclusions
The micropolar fluid over the shrinking surface has been considered. The system of governing equations has been transformed into the system of ODEs by using appropriate exponential similarity transformation. The system of ODEs is reduced to IVPs by employing the shooting method before solving IVPs by the Runge Kutta method. The pointwise conclusions of this study are given below:

2.
There are ranges of multiple solutions and no solutions that depend upon the suction parameter.

3.
According to stability analysis, the first solution is stable, which can be experimentally seen.

5.
The thickness of thermal boundary layer increases in the first and the second solutions as material parameter K is increased. 6.
Increments in the material parameter produce more couple stress. Funding: This research is also supported by the Universiti Utara Malaysia.