Mixed Convection Stagnation Point Flow of a Hybrid Nanoﬂuid Past a Permeable Flat Plate with Radiation Effect

: This article focuses on the stagnation point ﬂow of hybrid nanoﬂuid towards a ﬂat plate. The cases when the buoyancy forces and the ﬂow are in the opposite direction and the same direction are discussed. The effect of radiation and suction is also taken into account. The similarity trans-formations are used to convert the partial differential equations into nonlinear ordinary differential equations. These equations are computed numerically via the bvp4c function in MATLAB software. A comparison with the previously published articles is carried out, where an outstanding agreement is observed. The dual solutions exist in the case of opposing ﬂow ( λ < 0 ) and the suction parameter S > 0.6688. Meanwhile, only unique solutions exist in the case of assisting ﬂow ( λ > 0 ) . The existence of dual solutions leads to stability analysis. From the analysis, the ﬁrst solution is conﬁrmed as a stable solution. Furthermore, the heat transmission rate increases, while the skin friction coefﬁcient decreases as the radiation rate increases. An increase in the radiation rate from 0 (no radiation) to 1.0 increases the heat transmission rate by 5.01% for water, 4.96% for nanoﬂuid, and 4.80% for hybrid nanoﬂuid. Finally, it is worth mentioning that the present study yields new and original results. This study has also not been done by other researchers, indicating its novelty.


Introduction
Over the last few decades, the topic of nanofluid has attracted a vast number of researchers due to its contribution to industries. Note that the study of nanofluid flow was proposed by Choi and Eastman [1]. Nanofluid is a homogenous mixture of nanoparticles (i.e., metals, metal oxides, and carbon materials) and regular fluids (i.e., water, ethylene glycol, and oil). In contrast to regular fluids, nanofluids possess one-phase heat transfer coefficients and have larger thermal conductivity. Through the years, many investigations for development purposes have been conducted for enhancing heat transfer performance in certain flow problems. A new kind of fluid, namely hybrid nanofluid, is introduced by suspending two different nanoparticles in the regular fluid. By combining these nanoparticles, their chemical and physical properties will simultaneously combine and lead those properties in a homogeneous state [2,3]. The reviews on the preparation, thermophysical properties, and applications of hybrid nanofluid have been studied by several authors [4][5][6][7][8][9][10][11][12][13]. Recently, some published works have been reported for estimating the properties of both nanofluids and hybrid nanofluids. Esfe et al. [14] studied the thermal properties of Mg/O water nanofluid using artificial neural networks (ANNs). Later on, Fuxi et al. [15] investigated the thermal characteristics of water-EG/MWCNT-Al 2 O 3 hybrid nanofluid using a feed-forward neural network. The effect of nanoparticle volume fraction and temperature on dynamics viscosity of Al 2 O 3 -MWCNT (40:60)-Oil SAE50 hybrid nanofluid is examined by Qing et al. [16]. Esfe and Toghraie [17] applied an optimal feed-forward artificial neural network model and new empirical correlation to predict the viscosity of Al 2 O 3 -engine oil nanofluid. It is found that ANN estimated laboratory data more accurately than correlation output and ANN output.
Hybrid nanofluids have multitudinous applications in thermal energy transport, especially in manufacturing, naval structures, microfluidics, defense, transportation, acoustics, and biomedical [5,6,8]. Motivated by these important applications, numerous studies have been done to study the hybrid nanofluid in the boundary layer flow. For example, Devi and Devi [18] performed the hybrid nanofluid flow and heat transfer analysis over a stretching surface. They observed that the hybrid nanofluid had a higher heat transfer rate compared to nanofluid. Apart from that, Yousefi et al. [19] carried out an analytical study of the stagnation point flow near a cylindrical surface in a hybrid nanofluid. Their study showed that the thermal characteristics of the regular fluid are smaller compared to hybrid nanofluid and nanofluid with single nanoparticles. Moreover, Muhammad et al. [20] considered the stagnation point flow past a stretching sheet using cupric oxide (CuO) and carbon nanotubes (CNTs) with gasoline oil. Furthermore, continuous studies on hybrid nanofluid have been performed by many authors using various surfaces and physical situations [21][22][23][24][25][26].
The fluid motion near the stagnation region of a solid body is known as a stagnation point flow. The stagnation point flow can be defined as a point located at the object surface in the flow with zero local fluid velocity. In the stagnation region, the object leads the fluid to rest. When the fluid velocity is absent, the static pressure is at its highest magnitude at the stagnation point. The analysis of the stagnation point received much attention from researchers due to their important applications in some fields including nuclear reactor cooling systems, electronics cooling systems, aerospace, and so forth [27]. The stagnation point flow past a stationary semi-infinite plate was first discussed by Hiemenz [28], who applied the similarity variables to reduce the Navier-Stokes equations to nonlinear ordinary differential equations. Furthermore, consideration of mixed convection in stagnation point flow is one of the interesting topics on which to conduct research. In the presence of mixed convection, the buoyancy forces are caused by the temperature difference between the free stream and the surface enhancement. This leads to altering the flow as well as the thermal field. In this situation, the symmetrical behavior of the flow and thermal fields no longer exists concerning the stagnation line [29]. Other than that, this situation can also increase or decrease both local shear stress and heat transfer compared to forced convection flow.
Mixed convection flow near a stagnation point finds its applications in heat exchangers, electronic devices, solar and nuclear reactors, and atmospheric boundary layer flows. Comprehensive references on this topic can be found in the reported articles. Moreover, Tamim et al. [30] studied the MHD mixed convection stagnation point flow near a permeable surface in a nanofluid. In their study, the flow is subjected to prescribed surface temperature and external flow. Later on, Dinarvand et al. [31] examined the stagnation point flow of a nanofluid towards a vertical permeable stretching or shrinking plate. In addition, Rostami et al. [32] obtained two solutions for the mixed convective stagnation point flow of an aqueous silica-alumina hybrid nanofluid using the bvp4c solver. They found that the thermophysical characteristics of nanofluid and regular fluid can be enhanced by suspending two different nanoparticles. Finally, the study of combined free and forced convection flow in a stagnation region past a stretchable surface is considered by Seth et al. [33]. Thereafter, Zainal et al. [34] attempted to study the impact of convective boundary conditions and magnetic field on the flow towards a vertical flat surface in a hybrid nanofluid. The latest published work on mixed convection stagnation point flow is discussed by Ali et al. [35] under the influence of radiation and a magnetic field near a vertical stretching surface.
The boundary layer analysis with radiative heat flux is highly significant in the area of processes and space technology that requires high-temperature [36]. Radiation effects have numerous applications, for instance, in astrophysical flows, electricity generation, and solar power technologies. However, the addition of radiation in the thermal energy equation drives a highly nonlinear partial differential equation. Several works have been conducted to study the radiation effect in the boundary layer flow. For example, Turkyilmazoglu and Pop [37] performed an analytical study of a nanofluid flow adjacent to a vertical surface with the influence of radiation. The flow towards a permeable shrinking surface in a Casson fluid with inconstant surface temperature and radiation was examined by Bhattacharyya et al. [38]. Moreover, Soomro et al. [39] analyzed the influence of heat generation/absorption and radiation on the stagnation point flow near a stretching plate in a nanofluid. In addition, Jha and Samaila [40] numerically studied the thermal boundary layer flow on a flat surface with the convective boundary condition and thermal radiation effect. Furthermore, Anuar et al. [41] performed the unsteady micropolar hybrid nanofluid flow past a deformable sheet in the stagnation region. They noticed that the addition of radiation increased the heat transfer rate. Very recently, Jamaludin et al. [42] considered the mixed convection stagnation point flow in a cross fluid towards a permeable shrinking surface with radiation and suction effects.
Motivated by the above-mentioned articles, the novelty of the current work is to explore the existence of two solutions and the effect of radiation in mixed convection stagnation point flow near a permeable flat plate in Al 2 O 3 -Cu hybrid nanofluid. This problem follows the nanofluid correlations proposed by Devi and Devi [18]. It is worth mentioning that this problem has not been analyzed in any of the referenced state-of-the-art reviews before. To simplify the system of equations, the partial differential equations are converted into a non-dimensional form. The solutions for these equations are executed using the bvp4c package in the MATLAB program. Moreover, the use of this numerical method is inspired by the work of Rostami et al. [32], where several examinations are performed to figure out the accuracy of the current model. The current results are compared with several published articles for validation purposes. The physical reliability of the dual solutions gained is determined using stability analysis.

Description and Formulation of the Model
We considered a steady stagnation point flow of a hybrid nanofluid over a vertical flat plate with the effect of suction and radiation as shown in Figure 1. The x-axis is taken along the flat plate, while the y-axis is perpendicular to the surface of the plate. It is assumed that the ambient velocity of the flow is u e (x) = cx, where c is a constant with c > 0. The surface temperature is T s (x) = T ∞ + T * (x/L), where T ∞ is the ambient temperature, T * is the characteristic temperature and L is the characteristic length of the surface. The opposing flow (T * < 0) occurs when the upper part of the surface is cooled while the lower part of the surface is heated. In contrast, the assisting flow (T * > 0) occurs when the upper part of the surface is heated while the lower part of the surface is cooled. In this case, with the presence of buoyancy forces, the flow near the heated surface tends to move upward while the flow near the cooled surface tends to move downward. Hence, this behavior acts to assist the flow field [43]. Under the aforementioned assumptions, the governing equations are [30]: subject to [32] where (u, v) are the velocity fields of x and y directions. Respectively, g is the gravitational acceleration and T is the fluid temperature. The uniform mass flux is denoted by v 0 , where v 0 > 0 for injection and v 0 < 0 for suction. The radiative heat flux under Rosseland approximation is q r , which is given by [44]: In the above equation, σ * is the Stefan-Boltzman constant, while k * is the mean absorption coefficient. It is assumed that the temperature differences within the hybrid nanofluid are T 4 , which can be expressed as a linear function of temperature. Expanding the term T 4 in a Taylor series about T ∞ and neglecting higher-order terms, T 4 can be written as [44]: Thus, Equation (3) can be expressed as follows: The applied models for physical characteristics of hybrid nanofluid and nanofluid are given in Table 1. In the table, the subscript hn f , f , n f and s refers to hybrid nanofluid, fluid, nanofluid and solid nanoparticle, where s1 and s2 are the first and second nanoparticles, respectively. The nanoparticle volume fraction parameter is denoted as ϕ, where ϕ 1 is the nanoparticle volume fraction for alumina, and ϕ 2 is the nanoparticle volume fraction for copper. The thermophysical characteristics of copper (Cu), alumina (Al 2 O 3 ), and the regular fluid (water) are tabulated in Table 2. Table 1. Applied models for physical characteristics of hybrid nanofluid and nanofluid (see Devi and Devi [18]).

Hybrid Nanofluid
Dynamic viscosity, µ  Table 2. Physical properties of nanoparticles and the regular fluid [45]. The following similarity variables are introduced to convert our model into dimensionless form [30].

Hybrid Nanofluid
The uniform mass flux for the permeable surface is given by [32] v where S is the mass flux variable with S < 0 for injection and S > 0 for suction. Upon substituting Equation (8), Equation (1) is satisfied and Equations (2) and (7) become: where The transformed boundary condition can be written as: where prime is the differentiation with respect to ξ, Pr is the usual Prandtl number, Rd is the radiation parameter, and λ is the mixed convection parameter. The mentioned variables can be expressed as follows: Here, Re x = u e (x)x/ν f is the local Reynolds number and Gr is the local Grashof number. It is worth mentioning that, when ϕ 1 = ϕ 2 = 0 (regular fluid), Rd = 0 (radiation effects are negligible) and S = 0 (impermeable plate), Equations (10)- (12) reduced to those studied by Ramachandran [29]. The physical quantities of the skin friction coefficient, C f and the local Nusselt number, Nu x , are defined as follows [32]: where the shear stress of the plate, τ s , and the heat flux of the plate, q s , are given by [39]: Using Equations (8) and (15), Equation (14) then becomes:

Stability Analysis
Stability analysis is an analysis for identifying the stability of the solutions obtained. The similarity solutions emerge from Equations (10) and (11) subject to the boundary conditions (12). Note that several studies related to stability analysis have been done by other researchers [32,[46][47][48][49]. The procedures used in this analysis are referred to from the work of Merkin [50], Weidman et al. [51], and Harris et al. [52]. Following these three papers, the unsteady problem is initially considered as follows: Afterward, the dimensionless time variable, τ = ct, is introduced. Thus, the following new similarity variables are obtained: Considering terms from Equations (13) and (19), Equations (17) and (18) are converted into the following: 1 subject to the boundary conditions According to Weidman et al. [51], the similarity solutions are determined by substituting Equation (23) into Equations (20)- (22) with f = f 0 (ξ) and θ = θ 0 (ξ).
The linearized eigenvalue equations associated with the current problem are: together with the boundary conditions As discussed by Weidman et al. [51], τ is assumed to be zero indicating an initial decline or rise of the solution (23). Besides, the terms F(ξ, τ) and G(ξ, τ) can be expressed as F 0 (ξ) and G 0 (ξ), respectively. Harris et al. [52] stressed that it is necessary to relax one of the far field boundary conditions (26) whether on F 0 (ξ) or G 0 (ξ). In this work, we relaxed the condition F 0 (ξ) → 0 as ξ → ∞ and replaced it with a new condition F 0 (0) = 1.

Method of Solutions
Computation of results for Equations (10)-(12) is performed by using the bvp4c solver in MATLAB software. The solver is programmed with a collocation method with fourth-order accuracy (see Shampine et al. [53]). It is one of the effective techniques for solving boundary value problems (BVPs) for ordinary differential equations. The solver requires three kinds of information such as the equation to be solved, its related boundary conditions, and an initial guess for the solution. Mathematically, the bvp4c solver uses the finite difference method, where the output is obtained using an initial guess supplied at the starting mesh points. It is necessary to change the step size to obtain proper precisions. Since BVPs can have more than one solution, the program requires users to supply a guess for the solution desired. Before importing the model into the solver, Equations (10)-(12) must be converted into a system of first-order as follows [53]: and ya(1) = S, ya(2) = 0, ya(4) = 1, Here, a is the condition at the surface (ξ = 0) while b is the condition at the free stream (ξ = ∞). Next, the appropriate boundary layer thickness, initial guesses and inputs of the physical parameters (i.e., Rd, S, λ, Pr, ϕ 1 and ϕ 2 ) must be set to compute the desired solution. The numerical solution can be accepted when the boundary conditions ( f (ξ) → 1, θ(ξ) → 0 as ξ → ∞) are contented and there is an error or warning produced during the execution. The boundary conditions at ξ → ∞ are substituted by ξ = ξ ∞ . In this work, the appropriate value of ξ is chosen as ξ ∞ = 10, where ξ ∞ is situated outside the boundary layer thickness. The procedures are repeated to ensure the converged result secures a tolerance limit of 10 −6 [48]. The obtained results are then presented in tables and are plotted graphically.
Since we obtained two solutions, the stability analysis is performed to validate the stability of the solutions. In importing this analysis, the same procedures as above are applied using Equations (24)-(26) [53].
Based on these equations, the obtained result will represent a stable solution if the value of ζ is positive. Otherwise, the solution is unstable if the value of ζ is negative [51,52].

Results and Discussion
The effects of physical parameters on velocity, temperature, surface shear stress, local heat flux, skin friction coefficient, and local Nusselt number or heat transfer rate are discussed in this section. Three types of fluids namely regular fluid (ϕ 1 = ϕ 2 = 0), nanofluid (ϕ 1 = 0.02, ϕ 2 = 0) and hybrid nanofluid (ϕ 1 = ϕ 2 = 0.02) are also investigated. To validate our numerical procedure, the values of the shear stress f (0) and the local heat flux −θ (0) is compared with the previously published works by Tamim et al. [30] and Rostami et al. [32]. These comparisons can be found in Tables 3 and 4. As observed in the tables, the results are in excellent agreement with the above-mentioned references. In the work of Tamim et al. [30], they solved their problem by using the shooting method. The numerical method used by Rostami et al. [32] is similar to that used in the present study. Since the tables show a good agreement between the results, thus, we are confident in the use of the present method and model. Table 3. Comparison values of f (0) for multiple values of Prandtl number when Rd = S = ϕ 1 = ϕ 2 = 0 and λ = 1.

Tamim et al. [30]
Rostami et al. [ 2 and 3 demonstrate the surface shear stress and the local heat flux for both assisting (λ > 0) and opposing (λ < 0) flows and for several values of Rd when S = 1.0, ϕ 1 = ϕ 2 = 0.02 and Pr = 6.2. As depicted in the figures, larger Rd gives lower values of f (0) and −θ (0) compared to the case of no radiation (Rd = 0). In addition, dual solutions exist when λ < 0. Meanwhile, unique solutions are obtained when λ > 0. Here, λ c represents a turning point where the first and second solutions intersect. Beyond the turning point (λ < λ c ), the similarity solutions terminate. As the values of λ diminish from the positive (assisting flow) to negative (opposing flow) values, the shear stress and the local heat flux decrease for the first solution. In opposing flow, the buoyancy forces resist the fluid motion. Consequently, fluids in the boundary layer get retarded, acting as an adverse pressure gradient and reducing shear stress on the flat plate. Besides, assisting flow (heated plate) has a larger surface temperature compared to a fluid temperature, implying that the heat transfer takes place from the plate to the surrounding fluid.  Figures 4 and 5 illustrate the velocity and temperature distributions for various value of Rd when S = 1.0, λ = −8.0, ϕ 1 = ϕ 2 = 0.02 and Pr = 6.2. When Rd gets larger, the momentum boundary layer thickness thickens (see Figure 4). This leads to a reduction in the velocity field. The opposite trend is observed for the temperature field, which increases as Rd increases. Physically, radiation involves the emission or transmission of energy in the form of particles through a medium. Thus, increasing the parameter Rd tends to enhance the temperature of the flow and its thermal layer thickness.   Figures 6 and 7 present the values of f (0) and −θ (0) for three different types of fluids namely, water (ϕ 1 = ϕ 2 = 0), nanofluid (ϕ 1 = 0.02, ϕ 2 = 0) and hybrid nanofluid (ϕ 1 = ϕ 2 = 0.02) with suction when λ = −8.0, Rd = 1.0 and Pr = 6.2. The presence of a higher suction rate intensifies both values of f (0) and −θ (0). This is due to the fact that the high suction rate accelerates the random motion of base fluid particles in the flow. This criterion leads to higher shear stress. The suction effect also speeds up the transportation of heat from the plate to the fluid. This happens due to the permeability of the plate itself.
One can notice that water gives the largest shear stress and local heat flux compared to the others. Nanofluid and hybrid nanofluid contain more nanoparticles compared to water. Due to lack of space to collide, fluids with nanoparticles show a less favorable performance for shear stress formation. Besides, two solutions emerge when S > S c . Meanwhile, no solutions exist when S < S c .    x Nu x ) for some parameter Rd and for several types of fluids. The table shows that the skin friction coefficients decline with larger Rd. In the case of no radiation (Rd = 0), hybrid nanofluid takes the largest value of Re 1/2 x C f and Re −1/2 x Nu x . Meanwhile, in the presence of the radiation effect, water gives the largest skin friction coefficient. In addition to that, water exhibits the same behavior for the heat transfer rate when Rd is high enough (Rd = 1). As the radiation rate intensifies, more heat is transferred in the form of nanoparticles. Besides, increasing the radiation rate from 0 (no radiation) to 1.0 enhances the heat transfer rate by 5.01% for water, 4.96% for nanofluid, and 4.80% for hybrid nanofluid. The main contribution of the radiation effect in this study is to enhance the thermal performance of the flow. Based on the results obtained, we have achieved the requirement for the cooling processes. Table 6 presents the skin friction coefficient and the local Nusselt number for various values of S for several types of fluids. The values of S influence an increment in the skin friction coefficients. The presence of larger permeability of the plate accelerates the random motion of the suspended particles. As a result, they are increasing the friction on the boundary. Similarly, larger values of S tend to increase the heat transfer rate. Furthermore, water takes the largest skin friction coefficient when S is smaller (S = 1). Nevertheless, when S increases from S = 2 to S = 3, hybrid nanofluid shows a prominent outcome compared to others. Table 5. Values of Re 1/2 x C f and Re −1/2 x Nu x for Rd and several types of fluids when S = 1.0, λ = −8.0 and Pr = 6.2.

Physical Quantities
Rd Water Nanofluid Hybrid Nanofluid  Table 7 presents the minimum eigenvalues ζ of both solutions for some values of Rd and λ. The positive value of ζ signifies that the first solution is stable, while the result is contrary to the second solution. In addition, as λ approaches λ c , the values of ζ continuously decreases until it approaches 0. Based on the previous studies on stability analysis [51], there exists an initial decline of disturbances in the region of the stable flow.

Conclusions
The analysis of mixed convection stagnation point flow in a hybrid nanofluid past a permeable flat surface has been investigated under the influence of radiation. The resulting system of equations is solved numerically via the bvp4c function in MATLAB software. The outcomes revealed that two solutions were obtained in the opposing flow region. The solutions terminated when λ < λ c , where a unique solution was observed when λ > 0. The stability analysis verified that the first solution is stable, whereas the second solution is unstable. The skin friction coefficient decreased with a larger radiation effect. Meanwhile, the heat transfer rate intensified with the radiation effect. Increasing the radiation rate from 0 (no radiation) to 1.0 increased the heat transfer rate by 5.01% for water, 4.96% for nanofluid, and 4.80% for hybrid nanofluid. Besides, the skin friction coefficient and the heat transfer rate were enhanced with a higher suction rate. In the absence of radiation effect, hybrid nanofluid produced higher values of Re 1/2 x C f and Re −1/2 x Nu x . Meanwhile, in the presence of the radiation effect, water has the largest Re 1/2 x C f and Re −1/2 x Nu x . The skin friction coefficient and the heat transfer rate increased with the suction effect for three types of fluids. Moreover, water produced the largest skin friction coefficient with smaller S (S = 1.0). Meanwhile, hybrid nanofluid showed the prominent result for Re 1/2 x C f as S increased from S = 2 to S = 3. The significance of considering the radiation effect in this kind of flow is that it can reduce the skin friction coefficients for both assisting and opposing flows. Note that this criterion is good for the lubrication process. In addition, it can increase the heat transfer rate, which is suitable for certain applications that require high heat transfer rates (i.e., cooling processes).

Data Availability Statement:
The data presented in this study are available on request from the corresponding author. The data are not publicly available due to its large size.