Unsteady Stagnation Point Flow of Hybrid Nanofluid Past a Convectively Heated Stretching/Shrinking Sheet with Velocity Slip

Abstract

Unsteady stagnation point flow in hybrid nanofluid (Al₂O₃-Cu/H₂O) past a convectively heated stretching/shrinking sheet is examined. Apart from the conventional surface of the no-slip condition, the velocity slip condition is considered in this study. By incorporating verified similarity transformations, the differential equations together with their partial derivatives are changed into ordinary differential equations. Throughout the MATLAB operating system, the simplified mathematical model is clarified by employing the bvp4c procedure. The above-proposed approach is capable of producing non-uniqueness solutions when adequate initial assumptions are provided. The findings revealed that the skin friction coefficient intensifies in conjunction with the local Nusselt number by adding up the nanoparticles volume fraction. The occurrence of velocity slip at the boundary reduces the coefficient of skin friction; however, an upward trend is exemplified in the rate of heat transfer. The results also signified that, unlike the parameter of velocity slip, the increment in the unsteady parameter conclusively increases the coefficient of skin friction, and an upsurge attribution in the heat transfer rate is observed resulting from the increment of Biot number. The findings are evidenced to have dual solutions, which inevitably contribute to stability analysis, hence validating the feasibility of the first solution.

1. Introduction

The most common problem in boundary layer flow that had been treated so far is much focused on those for steady flows. Even though it is the steady flows that seem to have the utmost significance in real-world demands, some cases of time-varying in the boundary layer which is unsteady indicate an important role in several engineering problems. Some of the examples are start-up processes where the motions in rest are transits from a steady flow to another, and periodic motions of the working fluid [1]. The behavior of unsteady boundary layer flow describes an unusual pattern compared to the steady flow owing to the additional time-dependent terms in the governing equations, which exaggerated the separation of boundary layer and the fluid motion arrangement [2,3]. The thermal and mechanical properties of such an unsteady mechanism in the boundary layer approximation have been studied both analytically and numerically. Elbashbeshy and Bazid [4] presented the numerical investigation towards the unsteady stretching surface with heat transfer analysis. At the same time, Bhattacharya [5] managed to prove the existence of dual solutions in unsteady stagnation point flow towards a shrinking sheet by employing the shooting method approach coupled with a Runge–Kutta integration scheme. Bachok et al. [6] concluded that inclusion of the unsteadiness parameter offers a significant impact towards the boundary layer flow in nanofluid and Fan et al. [7] presented analytical solutions using the homotopy analysis method (HAM) and managed to advertise a highly precise analytical estimation which is in excellent agreement with the numerical results offered by the Keller box scheme. It is worth mentioning that a considerable amount of reviews on the unsteady stagnation point flow due to a stretching/shrinking surface have been accomplished by numerous researchers, including [8,9,10,11].

The stagnation point flow is one of the important topics in mechanics of fluid, in the way that stagnation point generally occurs in both engineering and science flow fields. The stagnation point flow could be identified in the extrusion process, polymer industry, and plane counter jet [12,13]. The ground-breaking research in this topic was first initiated by Hiemenz [14] who exposed an analytical explanation of two-dimensional stagnation point flow, and soon after, Homann [15] conducted a classical study of stagnation point in three-dimensional flow with regard to an axisymmetric case; whereas Howarth [16] tackled the problem of non-axisymmetric flow close to the stagnation area in three-dimensional analysis. Recently, Khashi’ie et al. [17,18], Fang and Wang [19], Waini et al. [20], and Zainal et al. [21] have scrutinized the stagnation point flow problems in diverse aspects with no-slip boundary conditions. Nevertheless, in numerous engineering occasions, the slip effect should be comprised, such as flow over lubricated or coated surfaces, rough or striated surfaces [22] and internal rare field gas flow [23]. Examples of industrial applications involving the slip boundary conditions are fluid flow on multiple interfaces, rare field fluid problems, and also the reacting flow in reactors [24,25]. Navier [26] and Maxwell [27] were the primary researchers who pioneered the study of linear slip boundary conditions, while Wang [28] has well reflected a comprehensive theoretical analysis considering the no-slip boundaries concentrating on the stagnation point flow. Rao and Rajagopal [29] have conducted an extensive evaluation and argument between the slip and no-slip condition, and Jusoh et al. [30] deliberated a modified nanofluid model towards a stretching/shrinking surface by considering a velocity slip parameter in three-dimensional flow. The study revealed that an increase in the velocity slip magnitude contributed to the intensification of skin friction coefficients.

In certain cases, the velocity slip or the non-adherence of the fluid to a solid boundary phenomenon was witnessed, for instance, in the micro-scale devices [31]. Fluids promoting slip are critical in technical applications such as polishing artificial heart valves and internal cavities [32]. The flow behavior and the shear stress in the fluid are rather distinctive with a slip at the wall boundary compared to those with no-slip condition. Besides that, the velocity slip effect does influence the heat transfer rate and was confirmed by Mukhopadhyay [33], who had investigated the slip impact of the unsteady mixed convective flow towards a porous stretching surface with heat transfer. The analysis found that the heat transfer rate declines with the velocity slip parameter, while it upsurges with the unsteadiness parameter. Mahapatra and Nandy [34] conducted a numerical study of the unsteady stagnation point flow past a shrinking sheet and heat transfer with the presence of slip effects in a viscous fluid. The results conveyed that with the increase of the velocity slip and unsteady parameter, the heat transfer rate is reported to escalate. Meanwhile, in nanofluid flow, Majumder et al. [35] specified that exertion of the partial velocity slip against the sheet surface is common. By relying on the finding in [35], Noghrehabadi et al. [36] then examined the impact of partial velocity slip on the nanofluid boundary layer flow and heat transfer past a stretching sheet. The work in [36] reported that an increment in the velocity slip effect decreases the momentum boundary layer thickness. Van Gorder et al. [37], who examined the nanofluid boundary layer flow over a stretching surface, also conveyed a similar result as [36] and further explained that no-slip condition is not applicable for fluid flows at nanoscales. Besides that, Dinarvand and Rostami [38] studied the rotating nanofluid flow and heat transfer with the presence of internal heating, velocity slip, and different shapes of nanoparticles. They showed that an increment in the velocity slip effect reduces the skin friction coefficient significantly. Researchers also tend to analyze the effect of velocity slip in the unsteady nanofluid flow as unsteady flow problems are more relatable to real-world applications. For instance, Seth et al. [39] studied the unsteady hydromagnetic nanofluid flow and heat transfer past a non-linearly stretching surface with Navier’s velocity slip and presented the analysis of entropy generation. Other valuable references regarding unsteady nanofluid flow with the velocity slip effect can be found in [40,41,42].

Ever since the outstanding inventions achieved by Choi and Eastman [43], who originated the brilliant idea of demonstrating the nanoparticle suspension in a base fluid and came out with the nanofluid term, a better type of working fluid is still being pursued. Acknowledging the sufficient improvement in the thermal conductivity of the conventional fluid is crucial, an advanced nanofluid form known as hybrid nanofluid is introduced, which intends to have highly developed heat conductivity. This modern type of fluid agent has fascinated numerous researchers owing to its reputation in the emergence and improvement of thermal characteristics in realistic applications, including micro-channel, heat pipes, heat exchangers, air conditioning systems, and mini-channel heat sink [44,45]. Gupta et al. [46] and Xian et al. [47] have reviewed the preparation method of hybrid nanoparticles along with the stabilization and its significance in industrial sectors. One of the critical elements in establishing a sustainable hybrid nanofluid suspension is selecting an appropriate combination of nanoparticles. The most widely used nanoparticles for the formation of hybrid nanofluid suspension are carbon materials (graphite, MWCNTs, CNTs), metals (Cu, Ag), metal oxides (Al₂O₃, CuO, Fe₂O₃), metal carbide, and a metal nitride. Madhesh and Kalaiselvam [48] conducted an experimental study to examine the features of hybrid nanofluid as a coolant agent, Tahat and Benim [49] examined the efficiency of hybrid nanofluid on flat plate solar collector, and they verified that the viscosity, thermal conductivity, and density of the working fluid had increased together with the concentration of Al₂O₃/CuO concentration, thus enhancing the solar collector proficiency. Some early research on hybrid nanofluid that employed the numerical method was done by Labib et al. [50], who investigated the impact of base fluids and hybrid nanofluid using a two-phase mixture model in forced convective heat transfer. Moghadassi et al. [51] revealed that the heat transfer performance is enhanced by adding the nanoparticles of Al₂O₃-Cu hybrid nanofluid while creating a small pressure drop in the system regime, Devi and Devi [52] focused on the mathematical inspection towards a stretching sheet. In contrast, the evaluation of heat transfers in the natural convection of Al₂O₃/water nanofluid and Al₂O₃-Cu/water hybrid nanofluid with a discrete heat source was explored by Takabi and Salehi [53]. Additional details on this topic are well described in the literatures [54,55,56].

To the best of the authors’ knowledge, the existing literature does not consider the unsteady stagnation point flow of hybrid nanofluids with the presence of velocity slip parameter and stability analysis in their models. Thus, the addressed issues above have inspired the authors to perform a numerical study in unsteady stagnation point flow towards a convectively heated stretching/shrinking sheet in alumina–copper/water (Al₂O₃-Cu/H₂O) with the impact of velocity slip on heat transfer. The hybrid nanofluid is recognized by dispersing Al₂O₃ nanoparticles into H₂O, followed by Cu with different volume fractions and the thermophysical properties of the hybrid nanofluid are adopted from Ghalambaz et al. [57] and Takabi and Salehi [53], which were based on the feasible physical assumptions and are in agreement with the conservation of mass and energy. The present work also utilized the bvp4c approach, which can be accessed in the MATLAB programming system towards solving the formulated problem. The existence of more than one solution is predictable; thus, an analysis of solution stability is completed to confirm the steadiness of the solutions which has an actual physical interpretation. The explanation of the results and the convergence of the obtained solutions are deliberated on in detail. Particular cases of current findings are evaluated in accordance with those of Mahapatra and Nandy [34] and Wang [58]. Furthermore, the consensus between previous and current findings is outstanding, and the agreement is excellent.

2. Mathematical Model

The unsteady two-dimensional stagnation-point flow of a hybrid Al₂O₃-Cu/H₂O nanofluid over a convectively heated stretching/shrinking sheet with the influence of velocity slip is considered in this research work, as illustrated in Figure 1 (see Dzulkifli et al. [59]). The stretching/shrinking velocity is denoted by $u_w\left(x,t\right)=bx/\left(1-ct\right),$ where $b$ denotes a constant corresponds to stretching $\left(b>0\right)$ and shrinking $\left(b<0\right)$ cases while $c$ signifies the unsteadiness problem and $u_e\left(x,t\right)=ax/\left(1-ct\right)$ is the velocity of the free stream where $a>0$ represents the strength of the stagnation flow. The ambient temperature and the reference temperature are $T_{\infty }$ and $T_0$, respectively. Now, we let the bottom of the sheet be heated by convection from a hot fluid at a specific temperature $T_f\left(x,t\right)=T_{\infty }+T_0\frac{a{x}^2}{2\nu }{\left(1-ct\right)}^{-3/2}$ which supplies a coefficient of heat transfer, expressed by $h_f.$ From all of the assumptions above; the governing boundary layer equations can be acknowledged as [34].

(1)$\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}=0,$

(2)$\frac{\partial u}{\partial t}+u\frac{\partial u}{\partial x}+v\frac{\partial u}{\partial y}=\frac{\partial u_e}{\partial t}+u_e\frac{\partial u_e}{\partial x}+\frac{{\mu }_{hnf}}{{\rho }_{hnf}}\frac{{\partial }^{2}u}{\partial y^2},$

(3)$\frac{\partial T}{\partial t}+u\frac{\partial T}{\partial x}+v\frac{\partial T}{\partial y}=\frac{k_{hnf}}{{\left(\rho C_p\right)}_{hnf}}\frac{{\partial }^{2}T}{\partial y^2},$

where $u$ denotes the component of velocity in $x-$ axis, $v$ is the velocity component in $y-$ axis, ${\mu }_{hnf}$ is the Al₂O₃-Cu/H₂O dynamic viscosity, ${\rho }_{hnf}$ the density of Al₂O₃-Cu/H₂O, $T$ is the Al₂O₃-Cu/H₂O temperature, $k_{hnf}$ is the thermal/heat conductivity of Al₂O₃-Cu/H₂O and ${\left(\rho C_p\right)}_{hnf}$ is the Al₂O₃-Cu/H₂O heat capacity. The boundary conditions, together with the partial slip for velocity, are set to

(4)$\begin{array}{l}u=u_w\left(x,t\right)+H_1\nu \frac{\partial u}{\partial y},v=0,-k_{hnf}\frac{\partial T}{\partial y}=h_f+\left(T_f-T\right)\mathrm{at}y=0,\\ u\to u_e\left(x,t\right),T\to T_{\infty }\mathrm{as}y\to \infty ,\end{array}$

where $H_1=H{\left(1-ct\right)}^{1/2}$ is the velocity slip factor, in which $H$ refers to the initial value of the velocity slip factor. The copper (Cu) thermophysical properties, along with aluminum oxide (Al₂O₃) and water (H₂O) nanoparticles, are provided in Table 1, as demonstrated by [60]. In the meantime, Table 2 issued the thermophysical properties hybrid nanofluid as established by [53,57]. The nanoparticles solid volume fraction is represented by $\varphi$, ${\rho }_f$ indicates the H₂O density, and ${\rho }_s$ is the density of the hybrid nanoparticle, $C_p$ is the constant pressure of heat capacity, while $k_f$ denotes the thermal conductivity of H₂O and $k_s$ is the hybrid nanoparticles thermal conductivity.

In order to express the governing Equations (1)–(3) concerning the boundary conditions (4) in a much simpler form, the subsequent similarity transformations are presented [34]

(5)$\psi ={\left(\frac{a\nu }{1-ct}\right)}^{1/2}xf\left(\eta \right),\theta \left(\eta \right)=\frac{T-T_{\infty }}{T_f-T_{\infty }},\eta ={\left(\frac{a}{\nu \left(1-ct\right)}\right)}^{1/2}y,$

where $\psi$ is the stream function that can be specified as $u=\partial \psi /\partial y,v=-\partial \psi /\partial y$ and $\eta$ is the similarity variable. Thus, we attain

(6)$u=\frac{ax}{\left(1-ct\right)}{f}^{\prime }\left(\eta \right),v=-{\left(\frac{a\nu }{1-ct}\right)}^{1/2}f\left(\eta \right).$

In view of the above relations, by employing the similarity variables (5) and (6), Equations (2) and (3) reduce to the following set of nonlinear similarity differential equations

(7)$\frac{{\mu }_{hnf}/{\mu }_f}{{\rho }_{hnf}/{\rho }_f}{f}^{‴}+f{f}^{″}-{f^{\prime }}^2+1-\epsilon \left(f^{\prime }+\frac{1}{2}\eta f^{″}-1\right)=0,$

(8)$\frac{1}{\mathrm{Pr}}\frac{k_{hnf}/k_f}{{\left(\rho C_p\right)}_{hnf}/{\left(\rho C_p\right)}_f}{\theta }^{″}+f{\theta }^{\prime }-2f^{\prime }\theta +\frac{\epsilon }{2}\left(\eta {\theta }^{\prime }+3\theta \right)=0.$

Here, $\epsilon$ measures the unsteadiness parameter with $\epsilon =c/a$, Pr represents the Prandtl number where $\mathrm{Pr}={\nu }_f/{\alpha }_f.$ Next, the initial and boundary conditions (4) now transform into

(9)$\begin{array}{l}f\left(0\right)=0,f^{\prime }\left(0\right)=\lambda +\gamma f^{″}\left(0\right),-\frac{k_{hnf}}{k_f}{\theta }^{\prime }\left(0\right)=\mathrm{Bi}\left[1-\theta (0)\right],\\ f^{\prime }\left(\eta \right)\to 1,\theta \left(\eta \right)\to 0,\mathrm{while}\eta \to \infty .\end{array}$

From Equation (9), $\lambda$ symbolises as the ratio of velocity parameter, $\gamma$ and Bi are the dimensionless velocity slip parameter and Biot number, respectively, which are described as

(10)$\lambda =\frac{b}{a},\gamma =H{\left(a\nu \right)}^{1/2},\mathrm{Bi}=\frac{h_f}{k_f}\sqrt{\frac{\nu \left(1-ct\right)}{a}}.$

Next, we define the skin friction coefficient $\left(C_f\right)$ and the local Nusselt number $\left(N{u}_x\right)$ as

(11)$C_f=\frac{{\tau }_w}{{\rho }_f{u}_e{}^2},N{u}_x=\frac{x{q}_w}{k_f\left(T_f-T_{\infty }\right)}.$

The shear stress along the $x-$ direction is represented by ${\tau }_w$, while $q_w$ signifies the surface heat flux that accentuated by

(12)${\tau }_w={\mu }_{hnf}{\left(\frac{\partial u}{\partial y}\right)}_{y=0},q_w=-k_{hnf}{\left(\frac{\partial T}{\partial y}\right)}_{y=0}.$

By applying (5) and (12) into (11), we acquire

(13)$\sqrt{{\mathrm{Re}}_x}{C}_f=\frac{{\mu }_{hnf}}{{\mu }_f}{f}^{″}\left(0\right),\frac{1}{\sqrt{{\mathrm{Re}}_x}}{\mathrm{Nu}}_x=-\frac{k_{hnf}}{k_f}{\theta }^{\prime }\left(0\right),$

provided that ${\mathrm{Re}}_x=\frac{u_{e}x}{{\nu }_f}$ is the local Reynolds number in x- axis.

3. Analysis of Solution Stability

By obeying the efforts of Merkin [61] and Merill et al. [62] from their outstanding discoveries of stability analysis scheme, the unsteady equations need to be deliberated in order to ultimately identify the reliable and stable solution since we notice the appearance of non-uniqueness solutions in the boundary value problem (7)–(9). Now, in accordance with the unsteady-state problem, a new similarity conversion is proposed

(14)$\begin{array}{l}u=\frac{ax}{1-ct}\frac{\partial f}{\partial \eta }\left(\eta ,\tau \right),v=-{\left(\frac{a\nu }{1-ct}\right)}^{1/2}f\left(\eta ,\tau \right),\theta \left(\eta ,\tau \right)=\frac{T-T_{\infty }}{T_f-T_{\infty }},\\ \eta =\sqrt{\frac{a}{\nu \left(1-ct\right)}}y,\tau =\frac{a}{1-ct}t.\end{array}$

Employing the similarity variables of Equation (14) to Equations (7) and (8), we now obtain the following converted differential equations

(15)$\frac{{\mu }_{hnf}/{\mu }_f}{{\rho }_{hnf}/{\rho }_f}\frac{{\partial }^{3}f}{\partial {\eta }^3}+\left(f+\frac{\epsilon }{2}\eta \right)\frac{{\partial }^{2}f}{\partial {\eta }^2}-{\left(\frac{\partial f}{\partial \eta }\right)}^2-\epsilon \frac{\partial f}{\partial \eta }-\left(1+\epsilon \tau \right)\frac{{\partial }^{2}f}{\partial \eta \partial \tau }+\epsilon +1=0,$

(16)$\frac{1}{\mathrm{Pr}}\frac{k_{hnf}/k_f}{{\left(\rho C_p\right)}_{hnf}/{\left(\rho C_p\right)}_f}\frac{{\partial }^2\theta }{\partial {\eta }^2}+f\frac{\partial \theta }{\partial \eta }-2\theta \frac{\partial f}{\partial \eta }-\frac{\epsilon }{2}\eta \frac{\partial \theta }{\partial \eta }-\frac{\epsilon }{2}3\theta -\left(1+\epsilon \tau \right)\frac{\partial \theta }{\partial \tau }=0,$

with respect to

(17)$\begin{array}{l}f\left(0,\tau \right)=0,\frac{\partial f}{\partial \eta }\left(0,\tau \right)=\lambda +\gamma \frac{{\partial }^{2}f}{\partial {\eta }^2}\left(0,\tau \right),-\frac{k_{hnf}}{k_f}\frac{\partial \theta }{\partial \eta }\left(0,\tau \right)=\mathrm{Bi}\left[1-\theta (0,\tau )\right],\\ \frac{\partial f}{\partial \eta }\left(\eta ,\tau \right)\to 1,\theta \left(\eta ,\tau \right)\to 0,\mathrm{as}\eta \to \infty .\end{array}$

In accordance with Weidman et al. [63], to test the stability of the steady flow $f\left(\eta \right)=f_0\left(\eta \right)$ and $\theta \left(\eta \right)={\theta }_0\left(\eta \right)$ which fulfil the boundary value problem and boundary conditions (refer to (7)–(9)), we write

(18)$f\left(\eta ,\tau \right)=f_0\left(\eta \right)+e^{-\omega \tau }F\left(\eta \right),\theta \left(\eta ,\tau \right)={\theta }_0\left(\eta \right)+e^{-\omega \tau }G\left(\eta \right),$

by which $\omega$ is the eigenvalue of unidentified variables, while functions $F\left(\eta \right)$ and $G\left(\eta \right)$ are relatively small to $f_0\left(\eta \right)$ and ${\theta }_0\left(\eta \right)$. The eigenvalue problems (15) and (16) result in an infinite group of eigenvalues ${\omega }_1<{\omega }_2<{\omega }_3......$ that detect an early decay when ${\omega }_1$ is positive, while an early growth of disruptions is observed when ${\omega }_1$ is negative, which exposes the unstable flow. Substituting (18) into (15)–(17), we develop

(19)$\frac{{\mu }_{hnf}/{\mu }_f}{{\rho }_{hnf}/{\rho }_f}\frac{{\partial }^{3}F}{\partial {\eta }^3}+\left(f_0+\frac{\epsilon }{2}\eta \right)\frac{{\partial }^{2}F}{\partial {\eta }^2}+F\frac{{\partial }^2{f}_0}{\partial {\eta }^2}-2\frac{\partial f_0}{\partial \eta }\frac{\partial F}{\partial \eta }+\left(\omega -\epsilon \right)\frac{\partial F}{\partial \eta }=0,$

(20)$\frac{1}{\mathrm{Pr}}\frac{k_{hnf}/k_f}{{\left(\rho C_p\right)}_{hnf}/{\left(\rho C_p\right)}_f}\frac{{\partial }^{2}G}{\partial {\eta }^2}+\left(f_0-\frac{\epsilon }{2}\eta \right)\frac{\partial G}{\partial \eta }-2\left({\theta }_0\frac{\partial F}{\partial \eta }+G\frac{\partial f_0}{\partial \eta }\right)+F\frac{\partial {\theta }_0}{\partial \eta }+\left(\omega -\frac{3}{2}\epsilon \right)G=0,$

and the boundary conditions are

(21)$\begin{array}{l}F\left(0,\tau \right)=0,\frac{\partial F}{\partial \eta }\left(0,\tau \right)-\gamma \frac{{\partial }^{2}F}{\partial {\eta }^2}\left(0,\tau \right)=0,-\frac{k_{hnf}}{k_f}\frac{\partial G}{\partial \eta }\left(0,\tau \right)-\mathrm{Bi}G\left(0,\tau \right)=0,\\ \frac{\partial F}{\partial \eta }\left(\eta ,\tau \right)\to 0,G\left(\eta ,\tau \right)\to 0,\mathrm{as}\eta \to \infty .\end{array}$

The heat transfer stability and steady-state flow solutions $f_0\left(\eta \right)$ and ${\theta }_0\left(\eta \right)$ was implemented via $\tau \to 0$, therefore $F=F_0\left(\eta \right)$ and $G=G_0\left(\eta \right)$ in (19)–(21). As a consequence, an early growth of Equation (18) is detected, and the subsequent generalized eigenvalue problem is recognized

(22)$\frac{{\mu }_{hnf}/{\mu }_f}{{\rho }_{hnf}/{\rho }_f}{F}_0^{‴}+\left(f_0+\frac{\epsilon }{2}\eta \right)F_0^{″}+F_0{f}_0^{″}-\left(2f_0^{\prime }-\omega +\epsilon \right)F_0^{\prime }=0,$

(23)$\frac{1}{\mathrm{Pr}}\frac{k_{hnf}/k_f}{{\left(\rho C_p\right)}_{hnf}/{\left(\rho C_p\right)}_f}{G}_0^{″}+\left(f_0-\frac{\epsilon }{2}\eta \right)G_0^{\prime }+F_0{\theta }_0^{\prime }-2\left({\theta }_0{F}_0^{\prime }+G_0{f}_0^{\prime }\right)+\left(\omega -\frac{3}{2}\epsilon \right)G_0=0,$

subject to

(24)$\begin{array}{l}{F}_0\left(0\right)=0,F_0^{\prime }\left(0\right)-\gamma F_0^{″}\left(0\right)=0,-\frac{k_{hnf}}{k_f}{G}_0^{\prime }\left(0\right)-\mathrm{Bi}{G}_0(0)=0,\\ F_0^{\prime }\left(\eta \right)\to 0,G_0\left(\eta \right)\to 0,\mathrm{as}\eta \to \infty .\end{array}$

The range of possible eigenvalues can be calculated by resting a boundary condition [64] of the present problem. In this study, we choose to repose ${F^{\prime }}_0\left(\eta \right)\to 0$, and the linear eigenvalue problems (22)–(24) are disclosed as ${F^{″}}_0\left(0\right)=1$ for a fixed value of ${\omega }_1$. It is worth mentioning that the values of ${\omega }_1$ are proficient in measuring the stability of the corresponding solutions $f_0\left(\eta \right)$ and ${\theta }_0\left(\eta \right)$.

4. Results and Discussion

The mathematical computations of this research work were achieved by employing the bvp4c function in the MATLAB programming system subject to the governing ordinary differential Equations (7) and (8) together with the boundary conditions (9). A bvp4c method is a notable tool for solving the boundary value problem that has been extensively established by various researchers to clarify the boundary value concern. In order to ensure the obtainment of desired solutions, early estimation of the primary mesh point and variations step size is crucial. Also, a reasonable assumption of the thickness of the boundary layer along with effective preliminary approximation is essential for defining the non-uniqueness solutions. The consistency of the results generated in the present study is evaluated with those in [34,58], as accessible in Table 3, which is in excellent agreement.

In this study, we recognized the hybrid Al₂O₃-Cu/H₂O nanofluid by dispersing the first nanoparticle ${\varphi }_1$ (alumina) into the base fluid (water) followed by the second nanoparticle ${\varphi }_2$ (copper) with various amounts of volume fractions. It is important to declare that the dispersion of discrete nanoparticle Al₂O₃/Cu is capable of developing Al₂O₃-Cu/H₂O and Cu-H₂O nanofluids. Apart from that, the Prandtl $\left(\mathrm{Pr}\right)$ number is set to be fixed at $\mathrm{Pr}=6.2$ corresponds to water as the reference-based fluid, while as for the hybrid Al₂O₃-Cu/H₂O nanofluid, the size of the nanoparticles is assumed to be standardized, and the thermophysical properties effect of nanoparticles agglomeration is ignored. The main purpose of the present study is to examine the influence of the control parameter such as the nanoparticles volume fraction $\left({\varphi }_1,{\varphi }_2\right)$, the unsteadiness parameter $\left(\epsilon \right)$, the velocity slip parameter $\left(\gamma \right)$, and the Biot number $\left(\mathrm{Bi}\right)$ towards the coefficients of skin friction variations $\left(f^{″}\left(0\right)\right)$ and the heat transfer rate $\left({\theta }^{\prime }\left(0\right)\right)$. The alumina and copper nanoparticles volume fraction in the present work is chosen within the range of $0.00\le {\varphi }_2\le 0.04$ which motivated by the experimental work done by Suresh et al. [65] who conducted the synthesis, characterization of ${\mathrm{Al}}_2{\mathrm{O}}_3-{\mathrm{Cu}/\mathrm{H}}_2\mathrm{O}$ nanocomposite powder for different volume concentrations 0.1%, 0.33%, 0.75%, 1%, and 2%. In their valuable study, the stability of the prepared nanofluids was determined by measuring the pH of nanofluids, and the nanofluid stability was found to diminish with increasing volume concentration. On another note, various values of the controlling parameter that has been used, are set within the following extent; $0.1\le \epsilon \le 0.2,$ $0.1\le \gamma \le 0.4,$ and $0.2\le \mathrm{Bi}\le 0.7$ to ensure the certainty of the obtained solutions. The present study also interested in witnessing the non-uniqueness solutions that appear in the ordinary governing differential of the ensuing problem, hence confirm the real and valid solutions through the stability analysis by utilizing the bvp4c approach in MATLAB operating system (MATLAB R2019b, MathWorks, Natick, MA, USA).

The generated results of Equations (7) and (8) perceive the non-uniqueness (dual) solutions together with the boundary conditions (9) to a specific scope of ${\lambda }_c$ where ${\lambda }_c$ manifests the non-uniqueness solutions meeting point, and this too is regarded as a critical point. There is a single solution at a critical point, and that is the distinct line. The flow separation happens to occur after the critical point, and the flow is no longer laminar and does not obey the principle of boundary layer theory. According to the numerical results attained in this research work, it is proven that a non-uniqueness solutions appear, namely first and second solutions, as depicted in Figure 2, Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9, Figure 10, Figure 11, Figure 12, Figure 13, Figure 14 and Figure 15. The existence of the non-uniqueness solutions contributes to the analysis of solution stability so that one may be able to verify the theoretically relevant solution. The first solution, however, is predicted to be reliable and fundamentally exist in practice. The smallest eigenvalues ${\omega }_1$ for some values of $\lambda$ when $\epsilon =0.1,{\varphi }_1={\varphi }_2=0.02,\mathrm{Bi}=0.2,$ and $\mathrm{P}\mathrm{r}=6.2$ with $\gamma =0.1,0.2,0.4,$ are tabulated in Table 4. The values of ${\omega }_1$ is noted approaching zero in the first and second solutions for $\lambda \to {\lambda }_c$ in selected cases of $\gamma .$ Hence, the authors may assure that the stability formulation and practices of the current problem are accurate and reliable. The flow is measured unstable if ${\omega }_1$ is negative because it implies an early development of disruptions that leads to flow separation. The smallest eigenvalue ${\omega }_1$ with positive value implies the flow is practically attainable and sustainable, which elucidates the solution stabilizing property to overwhelm the allowing disruptions. Also, it denotes an initial deterioration of disruptions that appear. The positive value of ${\omega }_1$ describes the stable mode of the flow as in the first solution. The significance of dual solutions that contribute to stability analysis is in revealing the other possibilities of flow behavior, which may be useful for future references in the extrusion process. For example, in the present study, we have presented that the second solutions are unstable solutions with negative eigenvalues. These unstable solutions with the negative eigenvalues infer the growth of the disturbance in the solutions, especially at some range of the stretching/shrinking parameter $\left(\lambda <-1.0\right).$ When the shrinking rate increases, it limits other external forces’ effect at the sheet and eventually showed the opposite behavior of the transport phenomena than the first solutions in the same range of $\lambda <-1.0,$ as the respective parameter varies. Moreover, any boundary value problem can generate more than one solution because of the nonlinearity in the boundary value problem; for example, see the mathematical model (7)–(9). Also, changes in the governing parameter values cause bifurcations in solutions that yield the existence of dual solutions [1]. Therefore, the mathematical model (7)–(9), which obeyed the boundary layer assumptions, managed to exhibit the variety in the fluid flow and heat transfer behavior through the uniqueness and existence of the solutions.

The influence of velocity slip on the skin friction coefficient and the local Nusselt number past a convectively heated stretching/shrinking sheet of hybrid Al₂O₃-Cu/H₂O nanofluid are displayed in Figure 2, Figure 3, Figure 4 and Figure 5. Figure 2 presents the coefficient of skin friction $\left(f^{″}\left(0\right)\right)$ towards the stretching/shrinking parameter $\lambda$, which revealed a decline of $f^{″}\left(0\right)$ with the occurrence of velocity slip effect $\left(\gamma =0.1,0.2,0.4\right)$ at the boundary. The same results are obtained in the previous literature, as reported by Dzulkifli et al. [59]. Based on a physical perspective, an improvement in the slip parameter reflects the fact that the vorticity produced by the stretching/shrinking velocity is gradually decreased, thus the vorticity stays restricted within the boundary layer for greater stretching/shrinking velocity with the same straining velocity of the stagnation flow, and subsequently, the steady solution is achievable for some broad values of $\lambda$, as stated by Mahapatra and Nandy [34]. Figure 2 emphasizes that the solution for a certain value of $\gamma$ persist prior to a critical value $\lambda =\left({\lambda }_c<0\right)$ at which the boundary layer splits from the convectively heated stretching/shrinking sheet and the solution on the basis of the boundary layer approximations is not feasible. In short, no solution exists for ${\lambda }_c<-1.6851$. Moreover, the expansion of $\gamma$ results in the increment of $\left|{\lambda }_c\right|$ suggesting that the velocity slip parameter is efficiently influenced by raising the range of dual solutions. Thus, it is proven that the existence of a slip velocity impact may prolong the separation of the boundary layer. Figure 2 also highlights as the sheet is stretching at the rate of $\lambda =1,$ the value of $f^{″}\left(0\right)=0$ which describes no frictional drag exerted at the convectively heated stretching/shrinking sheet.

Figure 3 exposes an upward trend in $-{\theta }^{\prime }\left(0\right)$ when the velocity slip parameter arises on the convectively heated stretching/shrinking sheet, which is proportionate to the heat transfer rate. Those findings are in line with the remarkable work reviewed by Mahapatra and Nandy [34] and Khashi’ie et al. [66]. Evidently, in the case of the first solution, the rate of heat transfer presents an ascendant trend with an increase of $\gamma$ while the second solution permits a top-down direction when the velocity slip $\gamma$ enlarges. From the existing and current findings, the authors can conclude that the velocity slip contributes to the improvement of the heat transfer rate significantly, prior to this case study. However, the authors would also like to declare that such results may vary if different control parameters are taken into consideration. The distribution of velocity $\left(f^{\prime }\left(\eta \right)\right)$ and temperature $\left(\theta \left(\eta \right)\right)$ profiles over several values of $\gamma$ in the case of a convectively heated shrinking sheet is certified in Figure 4 and Figure 5. As demonstrated in Figure 4, it is confirmed that as $\gamma$ increases, $f^{\prime }\left(\eta \right)$ is subsequently decreased in the first solution and revealed an upward trend in the velocity profile in the second solution with the presence of velocity slip condition. As the slip occurs, the velocity flow near the sheet in no longer equal to stretching velocity. This is due to the slip condition where the pulling force of the stretching sheet is partly shifted to the fluid, resulting in the decrement of the velocity profiles in the first solution. In short, the slip condition reduces the momentum transfer from the sheet to the fluid. A contradict results are obtained in $\theta \left(\eta \right)$ as the velocity slip effect is exaggerated, as shown in Figure 5. The first solution in the temperature profile distribution increases as $\gamma$ intensifies, while it diminishes in the second solution. This may occur due to the appearance of slip boundary conditions on the wall, which triggers the hybrid nanofluid to retain the velocity on the walls, and such slip may prevent from total heat exchange of the hybrid nanofluid. It is therefore noticed when the slip coefficient is applied, the difference in temperature rate is increased. Additionally, the effect of convective heating progress that has been reflected in the present study, which controls the temperature of the stretching/shrinking surface also might contribute to this phenomenon. As a matter of fact, greater convection leads to increased surface temperatures, which permit the thermal impact to penetrate deeper within the hybrid Al₂O₃-Cu/H₂O nanofluid.

Figure 6 and Figure 7 expose the coefficient of skin friction $\left(f^{″}\left(0\right)\right)$ and the rate of heat transfer $\left(-{\theta }^{\prime }\left(0\right)\right)$ of a conventional nanofluid $\left({\varphi }_1=0.02,{\varphi }_2=0\right)$ and hybrid nanofluid $({\varphi }_1=0.02,{\varphi }_2=0.02,0.04)$ past a convectively heated stretching/shrinking sheet when ${\varphi }_2$ varies from 0.00 to 0.02 where $\epsilon =0.1,\gamma =0.2,\mathrm{Bi}=0.2,$ and $\mathrm{Pr}=6.2$. Figure 6 manifests that addition in ${\varphi }_2$ which indicates the transformation of the conventional Al₂O₃-H₂O fluid to the hybrid Al₂O₃-Cu/H₂O nanofluid, upsurges the values of $f^{″}\left(0\right)$ once the sheet is shrinking. The viscosity of hybrid Al₂O₃-Cu/H₂O nanofluid rises when ${\varphi }_2$ increases, which eventually improves the fluid velocity over the convectively heated shrinking sheet, as proven in Figure 8. The velocity profile in Figure 8 clarifies that the momentum boundary layer thickness was diminished in response to the rise of ${\varphi }_2,$ thereby raising the velocity of the fluid and boosting the gradient of velocity. In fact, the thinner momentum boundary layer continues to evolve the wall shear stress as well as the convectively heated shrinking sheet, leading to an enhancement of $f^{″}\left(0\right).$ As $f^{″}\left(0\right)$ increases, the result implies the increase of the frictional drag exerted on the convectively heated shrinking surface, which may delay the boundary layer flow separation. Besides, Figure 6 also highlights when the sheet is stretching at the rate of $\lambda =1$, the value of $f^{″}\left(0\right)=0$ which explains no frictional drag exerted at the sheet surface. Meanwhile, Figure 7 illustrates an increasing trend of the heat transfer characteristic or $-{\theta }^{\prime }\left(0\right)$ when the values of ${\varphi }_2$ past a convectively heated shrinking sheet, and this trend holds true to the first solution but is in dispute with the results in the second solution. In essence, as the conventional Al₂O₃-H₂O nanofluid becomes the hybrid Al₂O₃-Cu/H₂O nanofluid, the heat transfer efficiency improves. This finding upholds the assumption of the convective heat transfer system can be improved by optimizing the nanoparticle concentration when ${\varphi }_2$ increased. The results obtained in Figure 6 and Figure 7 are consistent with Waini et al. [20] and Zainal et al. [21], whereby adding the concentrations of hybrid nanoparticles may contribute to the improvement of the heat transfer rate, accordingly. The temperature profile in Figure 9 describes the temperature variations when the conventional Al₂O₃-H₂O nanofluid becomes the hybrid Al₂O₃-Cu/H₂O nanofluid in both first and second solutions. The incline in the temperature of hybrid nanofluid proliferates the thermal conductivity, which may be triggered by the extra energy dispersed through the increment of the nanoparticles volume fraction over the state of convectively heated stretching/shrinking sheet.

Figure 10, Figure 11, Figure 12 and Figure 13 show the impact of the unsteadiness parameter $\left(\epsilon \right)$ towards a convectively heated stretching/shrinking sheet when $\epsilon$ shifts from 0.1 to 0.2. The hybrid Al₂O₃-Cu/H₂O nanofluid characteristic with regard to the coefficient of skin friction $\left(f^{″}\left(0\right)\right)$ is depicted as in Figure 10. Figure 10 captures that when the sheet shrinks, the increment in $\epsilon$ conclusively increases the trend of $f^{″}\left(0\right)$ in the first solution. An increment in the unsteadiness parameter results in the reduction of the boundary layer thickness, as depicted in Figure 12 and consequently upsurge the velocity gradient on the convectively heated stretching/shrinking sheet, thus $f^{″}\left(0\right)$ improves. Furthermore, the existence of nanoparticle volume fraction in the working fluid (Al₂O₃-Cu/H₂O) might also trigger the increment of $f^{″}\left(0\right)$ owing to an uplift of the hybrid nanofluid viscosity. This result is aligned with the study done by Ismail et al. [67]. Figure 11 presents the unsteadiness parameter effect towards the rate of heat transfer $-{\theta }^{\prime }\left(0\right)$. According to the generated results, the heat transfer rate increases when $\epsilon$ increases when the convectively heated sheet is shrinking. The forced convective heat transfer is indeed proportional to the effectiveness of nanofluid heat conductivity, hence raising the reduced local Nusselt number notably. The dimensionless velocity profiles $f^{\prime }\left(\eta \right)$ with a different value of $\epsilon$ are depicted in Figure 12, where the presence of dual velocity profiles is also observed. As illustrated in Figure 12, the first solution increases proportionally to the increment of $\epsilon$ values while the second solution displayed contradict results of the first solution, possibly because of the enhancement in the unsteadiness of the flow. Meanwhile, the same trend of the solution in Figure 12 also reflected the graph of temperature distribution $\theta \left(\eta \right)$ in the convectively heated stretching/shrinking sheet with the existence of unsteadiness parameter, as portrayed in Figure 13. Apart from that, it is proven that the second solution in both profiles, i.e., velocity and temperature distribution showed larger boundary layer thickness in each solution of the unsteady cases than those of the first solution.

Figure 14 and Figure 15 depict the variants of heat transfer rate $-{\theta }^{\prime }\left(0\right)$ and temperature distribution profile $\theta \left(\eta \right)$ with a different value of the Biot number $\left(\mathrm{Bi}=0.2,0.5,0.7\right)$ towards the convectively heated stretching/shrinking sheet. In Figure 14, the heat transfer rate shows an upsurge trend in the dual solutions along with the augmentation of Bi values. The Biot number signifies the conduction resistance ratio within the sheet to convection resistance at the sheet. The critical values of the different usage of Biot numbers suggest no substantial impact on the magnitude of $-{\theta }^{\prime }(0)$, as highlighted by Jusoh et al. [68]. An increment in the Biot number related to the improvement of convective heating is observed to decrease the fluid temperature proficiently in the first and second solutions, as displayed in Figure 15. This result is in contrast with the idea of a large Biot number representing larger internal thermal resistance of the sheet compare to the boundary layer thermal resistance since the temperature distribution profile spike as the value of Bi increases. However, note that the Biot number is specifically correlated to the coefficient of heat transfer $h_f$; therefore, it is conversely related to the thermal resistance of the current problem. Consequently, the heat resistance decreases as the Biot number increases, thereby increasing the heat transfer rate at the stretching/shrinking sheet and decreasing the temperature distribution (see Figure 15).

5. Conclusions

An analysis of the unsteady stagnation point flow of hybrid nanofluid over a convectively heated stretching/shrinking sheet incorporating the velocity slip impact on heat transfer was verified in this study. The results were generated by employing the bvp4c features in the MATLAB programming platform. The effect of diverse controlling parameters—namely, the nanoparticle volume fraction, the velocity slip, the Biot number, and the unsteadiness parameter—were reviewed. Our discoveries happen to determine that the presence of non-uniqueness solutions (first and second solutions) is verifiable within the hybrid Al₂O₃-Cu/H₂O nanofluid for a specific range of control parameters, and the stability analysis authorizes the reliability of our first solution. The nanoparticle volume fraction increment improved both the skin friction coefficient and local Nusselt number in the hybrid nanofluid (Al₂O₃-Cu/H₂O). From this, it is proven that the heat transfer rate improves when the ordinary Al₂O₃-H₂O nanofluid becomes the hybrid nanofluid (Al₂O₃-Cu/H₂O) by expanding the nanoparticles concentration. The recent addition of the velocity slip parameter at the boundary had encouraged a reduction in the skin friction coefficient and velocity profile. However, it increased the rate of heat transfer significantly. The temperature profiles escalate as the magnitude of velocity slip upsurge because of such slip may prevent from exchange of total heat in the hybrid nanofluid. An increase in the unsteadiness parameter consequently raises the velocity gradient on the convectively heated stretching/shrinking sheet, thus improving the skin friction coefficient. Besides, an upsurge in the Biot number intensity boosts the heat transfer rate since the Biot number is directly associated with the heat transfer rate coefficient. Thus, it is conversely interrelated to the thermal resistance of the designated problem. Apart from that, the critical values of the different practices of the Biot number imply no significant outcome towards the magnitude of the heat transfer coefficient.

Author Contributions

Research design, N.A.Z., R.N., K.N., and I.P.; Formulation and methodology, N.A.Z.; Result analysis, N.A.Z.; Validation, R.N. and K.N.; Article preparation, N.A.Z.; Review and editing, N.A.Z., R.N., K.N., and I.P. All authors have read and agreed to the published version of the manuscript.

Funding

The present work is endorsed by the research award (DIP-2017-009) from UKM.

Acknowledgments

All authors value the productive feedbacks by the competent reviewers.

Conflicts of Interest

The authors declare no conflict of interest.

Nomenclature

Roman letters
$a,b,c$	constant $(-)$
Bi	Biot number $(-)$
$C_f$	skin friction coefficient $(-)$
$C_p$	specific heat at constant pressure $\left({\mathrm{Jkg}}^{-1}{\mathrm{K}}^{-1}\right)$
$H_1$	velocity slip factor $(-)$
$h_f$	heat transfer coefficient $\left({\mathrm{Wm}}^{-2}{\mathrm{K}}^{-2}\right)$
$f\left(\eta \right)$	dimensionless stream function $(-)$
$k$	thermal conductivity of the fluid $\left({\mathrm{Wm}}^{-1}{\mathrm{K}}^{-1}\right)$
$N{u}_x$	local Nusselt number $(-)$
$\left(p{C}_p\right)$	heat capacitance of the fluid $\left({\mathrm{JK}}^{-1}{\mathrm{m}}^{-3}\right)$
$\mathrm{Pr}$	Prandtl number $(-)$
${\mathrm{Re}}_x$	local Reynolds number in $x-$ axis $(-)$
$t$	time $\left(\mathrm{s}\right)$
$T$	fluid temperature $\left(\mathrm{K}\right)$
$T_0$	reference temperature $\left(\mathrm{K}\right)$
$T_{\infty }$	ambient temperature $\left(\mathrm{K}\right)$
$u,v$	velocities component in the $x-$ and $y-$ directions, respectively $\left({\mathrm{ms}}^{-1}\right)$
$u_e$	velocities of the free stream in $\left({\mathrm{ms}}^{-1}\right)$
$u_w$	velocities of the stretching/shrinking surface $\left({\mathrm{ms}}^{-1}\right)$
x, y	rtesian coordinates $\left(\mathrm{m}\right)$
Greek symbols
$\psi$	stream function $(-)$
$\eta$	similarity variable $(-)$
$\theta$	dimensionless temperature $(-)$
$\epsilon$	unsteadiness parameter $(-)$
$\lambda$	ratio of the velocity parameter $(-)$
$\gamma$	velocity slip parameter $(-)$
$\mu$	dynamic viscosity of the fluid $\left({\mathrm{kgm}}^{-1}{\mathrm{s}}^{-1}\right)$
$\nu$	kinematic viscosity of the fluid (${\mathrm{m}}^2{\mathrm{s}}^{-1}$)
$\rho$	density of the fluid $\left({\mathrm{kgm}}^{-3}\right)$
$\tau$	dimensionless time variable $(-)$
${\tau }_w$	wall shear stress $\left({\mathrm{kgm}}^{-1}{\mathrm{s}}^{-2}\right)$
${\varphi }_1$	nanoparticle volume fractions for Al2O3 (alumina) $(-)$
${\varphi }_2$	nanoparticle volume fractions for Cu (copper) $(-)$
$\omega$	eigenvalue $(-)$
${\omega }_1$	smallest eigenvalue $(-)$
Subscripts
$f$	base fluid $(-)$
$nf$	nanofluid $(-)$
$hnf$	hybrid nanofluid $(-)$
$s1$	solid component for Al2O3 (alumina) $(-)$
$s2$	solid component for Cu (copper) $(-)$
Superscript
$\prime$	differentiation with respect to $\eta$ $(-)$

Figures and Tables

Figure 1. The schematic of problem flow (Dzulkifli et al. [59]).

Figure 2. Variants of f″(0) towards λ with γ=0.1,0.2,0.4.

Figure 3. Variants of −θ′(0) towards λ with γ=0.1,0.2,0.4.

Figure 4. Velocity profiles of f′(η) with λ=−1.3 (shrinking case).

Figure 5. Temperature profiles of θ(η) with λ=−1.3 (shrinking case).

Figure 6. Variants of f″(0) towards λ with ϕ2=0.00,0.02,0.04.

Figure 7. Variants of −θ′(0) towards λ with ϕ2=0.00,0.02,0.04.

Figure 8. Velocity profiles of f′(η) with λ=−1.4 (shrinking case).

Figure 9. Temperature profiles of θ(η) with λ=−1.4 (shrinking case).

Figure 10. Variants of f″(0) towards λ with ε=0.10,0.15,0.20.

Figure 11. Variants of −θ′(0) towards λ with ε=0.10,0.15,0.20.

Figure 12. Velocity profiles of f′(η) with ε=0.10,0.15,0.20.

Figure 13. Temperature profiles of θ(η) with ε=0.10,0.15,0.20.

Figure 14. Variants of −θ′(0) towards λ with Bi=0.2,0.5,0.7.

Figure 15. Temperature profiles of θ(η) with Bi=0.2,0.5,0.7.

Table 1

Cu thermophysical properties along with Al₂O₃ and H₂O (Oztop and Abu Nada [60]).

Properties	$\mathit{k}\left(\mathbf{W}/\mathbf{mK}\right)$	$\mathit{\rho }\left({\mathbf{kg}/\mathbf{m}}^3\right)$	${\mathit{C}}_{\mathit{p}}\left(\mathbf{J}/\mathbf{kgK}\right)$	$\mathit{\beta }\times {10}^{-5}\left(\mathbf{mK}\right)$
Cu	400	8933	385	1.67
Al2O3	40	3970	765	0.85
H2O	0.613	997.1	4179	21

Table 2

Hybrid Al₂O₃-Cu/H₂O nanofluids thermophysical properties (Takabi and Salehi [53], Ghalambaz et al. [57]).

Properties	Hybrid Nanofluid
Dynamic viscosity	${\mu }_{hnf}=\frac{1}{{\left(1-{\varphi }_{hnf}\right)}^{2.5}}$
Density	${\rho }_{hnf}=\left(1-{\varphi }_{hnf}\right){\rho }_f+{\varphi }_1{\rho }_{s1}+{\varphi }_2{\rho }_{s2}$
Thermal capacity	${\left(\rho C_p\right)}_{hnf}=\left(1-{\varphi }_{hnf}\right){\left(\rho C_p\right)}_f+{\varphi }_1{\left(\rho C_p\right)}_{s1}+{\varphi }_2{\left(\rho C_p\right)}_{s2}$
Thermal conductivity	$\frac{k_{hnf}}{k_f}=\left[\frac{\left(\frac{{\varphi }_1{k}_{s1}+{\varphi }_2{k}_{s2}}{{\varphi }_{hnf}}\right)+2k_f+2\left({\varphi }_1{k}_{s1}+{\varphi }_2{k}_{s2}\right)-2{\varphi }_{hnf}{k}_f}{\left(\frac{{\varphi }_1{k}_{s1}+{\varphi }_2{k}_{s2}}{{\varphi }_{hnf}}\right)+2k_f-\left({\varphi }_1{k}_{s1}+{\varphi }_2{k}_{s2}\right)+{\varphi }_{hnf}{k}_f}\right]$

Table 3

Evaluation values of $f^{″}\left(0\right)$ when $\epsilon =\gamma =\mathrm{Bi}=0,$ by certain values of $\lambda$.

$\mathit{\lambda }$	Present Result		Mahapatra and Nandy [34]		Wang [58]
	First Solution	Second Solution	First Solution	Second Solution	First Solution	Second Solution
−0.25	1.402241	-	1.402242	-	1.4022404	-
−0.50	1.495670	-	1.495672	-	1.4956704	-
−0.75	1.489298	-	1.489296	-	1.4893004	-
−1.00	1.328817	0.000000	1.328819	0.000000	1.3288204	0.000000
−1.10	1.186680	0.049229	1.186680	0.049229	-	-
−1.15	1.082231	0.116702	1.082232	0.116702	1.082230	0.116702
−1.20	0.932473	0.233650	0.932470	0.233648	-	-
−1.246	0.609826	0.529035	0.584374	0.554215	0.5543004	-

Table 4

Smallest eigenvalues ${\omega }_1$ for some values of $\lambda$ when $\gamma =0.1,0.2,0.4$.

$\mathit{\gamma }$	$\mathit{\lambda }$	${\mathit{\omega }}_1$First Solution	${\mathit{\omega }}_1$Second Solution
0.1	−1.30	0.7646	−0.7093
	−1.380	0.1015	−0.1398
	−1.3820	0.0095	−0.0500
0.2	−1.40	0.6870	−0.6316
	−1.460	0.2387	−0.2435
	−1.4690	0.0088	−0.0225
0.4	−1.68	0.1822	−0.1576
	−1.684	0.0912	−0.0696
	−1.6851	0.0276	−0.0067