An efficient numerical method for 2D elliptic singularly perturbed systems with different magnitude parameters in the diffusion and the convection terms, part Ⅱ

1.   Introduction

A significant value has been bestowed to the problems containing non-Newtonian fluids because it has a wide range of applications in industries and engineering fields. It is frequently experienced in real-life problems, which captivated remarkable attraction from scholars owing to their various significances in the fields like polymer melts, drilling mud, blood polymers, fruit juices, certain greases and oils, suspensions, etc. It is intricate to propose one model that contains all rheological properties of non-Newtonian fluids. Different models like Maxwell fluid ^[1], Burger's fluid ^[2] were suggested by scholars to overcome these difficulties. Williamson fluid model is one of the best models among these non-Newtonian fluids because it involves both maxima as well as minimum viscosities that provide better results for pseudo-plastic fluids. This model was proposed in 1929 by Williamson ^[3], which explains the viscous flow equations for pseudo-plastic fluids and also results verified experimentally. Nadeem et al. ^[4,5] analyzed 2D flow conveying Williamson fluid over linear as well as exponential stretching sheets and utilized HAM to obtain a series solution. Khan and Khan ^[6] hired the HAM procedure to find the result of Williamson fluid. Hayat et al. ^[7] inspected the impact of the magnetic field on the unsteady flow of Williamson fluid through a permeable plate and observed that the Williamson parameter decelerates the velocity gradient. The impact of chemical reaction and features of heat transport of Williamson fluid over a stretching surface with nanofluid was scrutinized by Krishnamurthy et al. ^[8]. Kumar et al. ^[9] studied 3D flow conveying Prandtl fluid through a heated Riga surface subject to a chemical reaction. Zaib et al. ^[10] explored the flow of time dependent function comprising nanofluid of Williamson fluid over a heated porous shrinking surface and found dual solutions for decelerated flow. Khan et al. ^[11] found multiple solutions of nanofluid flow comprising Williamson fluid over a moving thin needle with erratic radiation. Qureshi ^[12] inspected the entropy generation on a steady flow of Williamson fluid conveying the water-based TiO₂ nanoparticle via a slippery surface with radiation and magnetic effects. Recently, Raza et al. ^[13] inspected the impact of thermal radiation on the fluid flow and phenomena of heat transfer induced by Williamson fluid via a curved surface. They establish that the velocity distribution declines due to the Williamson parameter, whereas the temperature of fluid uplifts.

In recent times, hybrid nanofluids have been extensively utilized by several scholars and researchers because hybrid nanofluids have superior thermal conductivity. Hybrid nanofluid is generated by either scattering dissimilar nanomaterials in composite or blend form. Hybrid nanofluid is utilized in several engineering and industrial fields like naval structures, nuclear safety, pharmaceutical, cooling in modern electronic devices, heat transport, and drug reduction. Jana et al. ^[14] compared the enhancement of thermal conductivity of hybrid and single nano-additives. Similar research has been conducted by Jha and Ramaprabhu ^[15] through employing hybrid nanofluid as well as MWCNT and observed that thermal conductivity drastically increased compared to single nanofluid. The properties of HN data experimentally as a mechanism of coolant was scrutinized by Madhesh and Kalaiselvam ^[16], whereas Devi and Devi ^[17,18] conducted an analysis mathematically with water-based Cu-Al₂O₃ hybrid nanofluid and heat transport via a stretched sheet. In addition, the hybrid nanofluid has been extensively reviewed by several scholars ^{[19,20,21,22,23]}. Recently, Sheikholeslami et al. ^[24] inspected the impacts of radiation and Lorentz forces conveying hybrid nanofluid in a non-Darcy medium and utilized the CVFEM technique to obtain the solution. The steady flow containing hybrid carbon nanotubes over a stretched permeable cylinder with magnetic consequence is examined by Gholinia et al. ^[25]. Khan et al. ^[26] studied the mixed convection flow containing ethylene-based HN over a vertical shrinking or stretching cylinder and stability analysis was performed. Waini et al. ^[27] used the magnetic influence on a steady flow with features of heat transport comprising HN through a shrinking/stretching wedge and found dual solutions. Wakif et al. ^[28] explored the surface roughness and radiation effects of conveying water-based copper-oxide and alumina nanoparticles over an infinite horizontal enclosure. The stagnation point flows over a stretching/shrinking surface in a hybrid nanofluid with various physical situations were considered by Zainal et al. ^[29] and Waini et al. ^[30], where the existence of dual solutions was reported. Bhatti and Sara ^[31] discussed the motion of the peristaltic flow of a Carreau fluid comprising a hybrid nanofluid through magnetic effects. They observed that targeting drug delivery is an accurate and effective way for drug delivery. The features of hybrid nanofluid with heat transfer unruffled by ZnO-Al₂O₃ particles with timed-dependent magnetic field were inspected by Zhang et al. ^[32]. They presented double solutions for the time-dependent parameter. Zhang et al. ^[33] examined the dynamics of magneto nanofluid through an erratic permeable stretchable sheet with wall transpiration and Arrhenius chemical reaction. They found that the thermophoretic parameter elevates nanoparticle concentration as well as temperature. Recently, Bhatti et al. ^[34] investigated the magnetic flow of a non-Newtonian nanofluid through circular rotating plates saturated in a porous media with motile organisms and utilized the DTM-Pade technique to find the numerical results.

The exploration of cross-flow was instigated for a time after the pioneering investigations through the papers by Prandtl ^[35] and Blasius ^[36] comprising the laminar flow from a flat surface by tinny viscosity. Prandtl ^[37] is the first apparently to report the outcome for homogeneous pressure gradient flow over a yawed finite cylinder. Cooke and Hall ^[38] and Eichelbrenner ^[39] respectively reviewed numerous theoretical as well as numerical techniques development for working out generalized 3D boundary-layer flows. Weidman ^[40] scrutinized boundary-layer through cross flows compelled via transverse motions of the plate. Recently, Roşca et al. ^[41] extended the work of Weidman by considering rotated stagnation point flow conveying hybrid nanofluid through a porous shrinking or stretching surface.

Consideration of thermal radiation adds a new dimension to boundary layer flow and heat transfer. The effect of radiative heat transfer on dissimilar flows is critical in high-temperature processes and space technology. Thermal radiation can have a substantial impact on heat transfer and distribution of temperature in the boundary layer flow of a participating fluid at high temperatures. Significant efforts have been made in recent years to achieve better control over the rate of cooling. Viskanta and Grosh ^[42] studied the thermal radiation effects on the flow and heat transfer past a wedge in emitting and absorbing media. The impact of thermal radiation on the heat transfer across the stretchable sheet was investigated by Elbashbeshy ^[43]. El-Aziz ^[44] and Hayat et al. ^[45] considered the consequence of radiation (RN) on the time-dependent flow and heat transfer past a stretched sheet. Zaib and Shafie ^[46] inspected the impression of RN on the time-dependent magneto flow past a stretched sheet with Joule heating, thermal stratification, and Hall current. Recently, Khan et al. ^[47] reviewed the RN effect on the blood flow of Williamson fluid induced by nanofluid past porous curved surface and present dual solutions.

Researchers have been studying the model of BLF and heat transfer over a shrinking/stretching surface (SN/ST) for several eras of time. This is because of its numerous and significant applications in technological and industrial queries such as wire drawing, aerodynamic extrusion of plastic sheets, hot rolling, metal spinning and so on. Sakiadis ^[48] pioneered the steady 2D viscous fluid flow past a moving solid continuous surface. Crane ^[49] became the first who investigate the steady 2D flow through a linearly stretchable sheet and discover the similarity outcome in analytical forms. Moreover, the emergence of an unconscious category of flow as a result of SN when studied the pattern of a liquid film on a time-dependent stretchable sheet was innovated by Wang ^[50]. This flow towards a SNS is effectively a reverse flow as argued by Goldstein ^[51]. Miklavčič and Wang ^[52] stated that transparent of mass is needed to maintain flow past a SNS. The flow stimulated by a SNS exhibits physical processes that differ significantly due to the onward elongating flow addressed by Fang et al. ^[53]. Since then, numerous scholars ^[54,55,56] recently inspected the stretching/shrinking flows with different aspects.

It is observed from current literature that the cross-flow induced by non-Newtonian fluid past a stretching/shrinking sheet is not yet investigated near a rotating stagnation point in the presence of shape factors. Thus, the novelty of the problem is to inspect the cross-flow of Williamson fluid conveying hybrid nanofluid over a shrinking/stretching sheet with radiation effect. The similarity multiple outcomes were found for a certain range of shrinking/stretching parameters. It is believed that the achieved numerical solutions employing the current suggested model for hybrid nanofluid set a novel scope for scholars in the analysis of heat transfer.

2.   Fluid model

The Cauchy stress tensor ${S}_{1}$ for the Williamson fluid is defined as

where ${\tau }_{A}$, ${\mu }_{0}$, ${\mu }_{\infty }$ signify as extra stress tensor, limiting viscosity at zero and infinite shear stresses, respectively. ${\varGamma }_{1} > 0$ represents the time, ${A}_{1}$ indicates the first Rivlin-Erickson tensor, and $\ddot{\gamma }$ is given as

Here, we discuss the case for which

We get

or

3.   Mathematical modeling of the flow problem

A two-dimensional incompressible, viscous, and laminar BL cross-flow and heat transfer of a Williamson hybrid nanofluid is considered over a porous shrinking/stretching surface. In existing scrutiny, the various shape factors, as well as thermal radiation are taken into account. The HN is made up of two changed nanoparticles (copper oxide CuO and titanium dioxide TiO₂) mixed together with a conventional base fluid (water H₂O). The properties of the hybrid nanofluid are taken to be in thermal equilibrium. The summarized picture of the flow configuration model along with the nanoparticles shape factor in the Cartesian coordinate system is shown in Figure 1, where the coordinate $x-$ axis is run along the surface of the sheet and the coordinate $y-$ axis is normal to it. The schematic diagram is such that the flow is protracted infinite in the coordinate $z-$ axis (span-wise) and thus, it is completely established. So we look for multiple outcomes with velocity profiles free of the coordinate $z-$ axis. As the surface of the sheet is stretched and contracted in the positive and negative $x-$ directions with a variable velocity such as ${u}_{w}\left(x\right) = {{\beta }_{A}}^{2/3}{\nu }_{Bf}^{1/3}{x}^{1/3}{\lambda }_{A}$, where ${\lambda }_{A}$ is the uniform stretching/shrinking parameter. Moreover, the surface of the sheet is stretching and shrinking owing to the positive and negative values of the parameter ${\lambda }_{A}$, respectively, while the zero value corresponds to the static sheet.

Additionally, supposed the non-uniform velocity of the free-stream or ambient from the sheet surface (inviscid flow) is represented by (Weidman et al. ^[57])

where ${\beta }_{A}$ signifies an arbitrary positive constant. Also, supposed that ${v}_{w}\left(x\right)$ corresponds to the transpiration or permeable velocity while ${v}_{w}\left(x\right) < 0$ stands for injection or blowing and ${v}_{w}\left(x\right) > 0$ stands for suction, respectively. The surface constant temperature and temperature of the far-field (inviscid flow) are denoted by ${T}_{w}$ and ${T}_{\infty }$, respectively. Under the aforementioned stated proviso, the requisite governing equations for the considered heat and fluid flow dynamics are the following as follow ^[11,12]:

For the current problem, the physical boundary conditions ^[57,58] are

Here, the components of velocity in the corresponding coordinates of $y$ and $x$ are characterized with $v$ and $u$, respectively, ${\varGamma }_{1}$ is the time constant and $T$ is the temperature of the fluid. Additionally, ${\left(\rho {c}_{p}\right)}_{Hbnf}$ signifies the heat capacity of hybrid nanofluid, ${k}_{Hbnf}$ signifies the thermal conductivity of the hybrid nanofluid, ${\rho }_{Hbnf}$ signifies the density of the hybrid nanofluid, and ${\mu }_{Hbnf}$ is the dynamic viscosity of the hybrid nanofluid.

In the aforementioned governing equations, the expression and constraints for the hybrid nanofluid are demarcated as (Nisar et al. ^[58] and Khan et al. ^[59]):

where ${c}_{p}$ the specific heat at constant pressure, ${k}_{Bf}$, ${\mu }_{Bf}$, ${\rho }_{Bf}$ and ${\sigma }_{Bf}$ stand for the respective thermal conductivity of the base fluid, the dynamic viscosity of the base fluid, the density of the base fluidand electrical conductivity of the base fluid. The subscript $Bf$, $\mathrm{C}\mathrm{u}\mathrm{O}$ and $\mathrm{T}\mathrm{i}{\mathrm{O}}_{2}$ correspond the respective quantities of the base fluid, copper oxide (CuO) nanoparticles and titanium dioxide (TiO₂) nanoparticles. Therefore, Table 1 encloses the data of the working pure fluid and two distinct nanomaterials (copper oxide (CuO) and titanium dioxide (TiO₂)). Moreover, the various shape factor coefficient and distinct shapes of nanoparticles are highlighted in Table 2.

In energy equation (3.4), the last term represents the thermal radiation, where ${q}_{rA}$ identifies the radiative heat flux and it is demarcated by exercising the Rosseland estimation ^[12], we have

where ${k}_{A}^{*}$ corresponds to the coefficient of mean absorption and ${\sigma }_{A}^{*}$ called the Stefan-Boltzmann constant. Now employing the Taylor series for the term of ${\left(T\right)}^{4}$ at a position ${T}_{\infty }$ and ignoring the higher order terms in approximation, which can take place the following final form as follow:

In the existing problem to straightforward further our model analysis by considering the subsequent non-dimensional variables (Weidman et al. ^[57]):

Now the primes denote derivatives with respect to the pseudo-similarity variable $\xi$ and the non-uniform velocity at the porous shrinking/stretching surface of the sheet is demarcated as:

Here, ${S}_{A}$ represents the factor of mass flux through ${S}_{A} > 0$ and ${S}_{A} < 0$ correspond to suction and blowing, respectively, whereas ${S}_{A} = 0$ used for an impermeable surface of the sheet.

Now utilizing the similarity transformations equation (3.9) into governing equations (3.2)–(3.4), where the continuity equation (3.2) is fully verified while the other requisite equations take place the following dimensionless form:

with the subjected BCs:

According to Weidman et al. ^[57], for $\alpha = 1$, Eq (6.1) is the same as Eq (3.11) when ${\varphi }_{Ti{O}_{2}} = {\varphi }_{CuO} = 0$ and $W{e}_{A} = 0$, while Eq (3.12) has not been implemented by ^[57]. Also, the involvement of dimensionless constraints in Eqs (3.11) and (3.12) are namely called the Weissenberg parameter, Prandtl number and radiation parameter. These parameters are mathematically expressed as

The physical quantities of engineering practical interest are the shear stress and rate of heat transfer which are demarcated as follows:

Utilizing Eq (3.9) into Eq (3.15), one obtains the following dimensionless form as

Finally, the term ${\mathit{Re}}_{x} = \frac{x{u}_{w}}{{\nu }_{Bf}}$ indicates the local Reynolds number.

4.   Numerical technique

This portion shows the solution procedure of considered flow and analysis of heat transfer suspended HN. After utilizing the self-similarity transformations (3.9), the model requisite equations are bounded in the form of highly nonlinear ordinary differential equations (3.11) and (3.12) along with the boundary conditions (3.13). These equations are solved numerically using a built-in function called bvp4c accessible in a MATLAB software (see Shampine et al. ^[61,62], Chu et al. ^[63]). It should be revealed that this bvp4c package is set up based on the finite difference scheme, which is additionally highlighted via the three-stage Lobatto IIIA formula. To execute the bvp4c technique, the higher third order and second order achievable non-dimensional form of ordinary differential equations (ODEs) is altered into a system of first-order ODEs by introducing new variables. By working out this procedure, let

Now substituting these variables in the achievable highly non-linear dimensionless form of ODE's along with BCs to get the transformed nonlinear ODE's in first order as follow:

with the initial conditions

For solving the system of Eq (4.2) along with the respective conditions (4.3), the method required initial early guesses at the mesh point. A continuous outcome is obtained using the polynomial such as the collective-type formula. The outcome is provided by a fourth-order precision set that is equivalently scattered over the spatial intervals where the function is integrated. According to many fruitful practices with the subject of boundary layer theory, the limitation at far away ($\xi \to \infty$) is substituted by the value $\xi = {\xi }_{\infty } = 8$ and the relative error tolerance is pre-demarcated as 10^-6. In addition, the error control and the mesh selection are based on the residual of the smooth outcome. The beginning mesh has 4 points equal discretized on the interval $\left[0, {\xi }_{\infty } = 8\right]$ and then the selection of mesh is adjusted automatically via the package of bvp4c. Currently, the problem has dual outcomes which mean that the package of bvp4c required two different guesses for branch of stable and unstable outcomes. The initial early guess for the branch of upper outcome is simple while the guess selection for the branch of lower outcome is quite difficult. According to Merkin et al. ^[64] and Weidman et al. ^[65], the branch of the upper outcome is physically reliable and stable while the branch of the lower outcome is unstable and not physically reliable because the outcome exists only for a particular space of the SN/ST parameter. Similar outcomes may arise in other circumstances, where the corresponding results may have more practical meaning (see Ridha ^[66] and Ishak et al. ^[67]). Finally, for the details and better interpretation of the given numerical scheme a flow chart (Figure 2) has also been added.

To substantiate the validity, accuracy and precision of the considered numerical scheme, the results of the shear stress of the given study for the limiting case of the shrinking parameter ${\lambda }_{A} < 0$ have been matched with the consequences of Waini et al. ^[68]. The numerical outcomes are shown in the posited Table 3, which displays an exceptional matching and harmony with the available outcomes that give us more assurance that the outcomes accessible in the existing paper are very precise.

5.   Results and discussion

This segment of the work demonstrates the two different solution branches for a certain range of the shrinking parameter ${\lambda }_{A}$. Results are bounded in the form of various pictures for the water-based copper oxide-titanium dioxide hybrid nanofluid from Figures 3–17 owing to the impact of various involved physical parameters like volume fractions of nanoparticles ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$, suction parameter ${S}_{A}$, Weissenberg number $W{e}_{A}$ and radiation parameter ${R}_{Ad}$ on velocity profile, temperature distribution, shear stress $f\text{'}\text{'}\left(\xi \right)$ in the whole boundary layer region, the skin friction and rate of heat transfer. For the computation, we have fixed the default values of the constraint parameter throughout the manuscript are the following such as ${\varphi }_{CuO} = 0.025$, ${\varphi }_{Ti{O}_{2}} = 0.035$, ${S}_{A} = 2.0$, $W{e}_{A} = 0.5$, $M = 3.0$ and ${R}_{Ad} = 5.0$. The value of the Prandtl number is taken to be 6.2. Following the whole computations, we observed that the impact of the shape factors on the gradients and velocity profiles are very diminishing. Therefore, we have drawn only one figure for the influence of various shape factors on temperature profiles while the rest of the figures were prepared for the case of the sphere $M = 3.0$.

Figure 3 illustrates the effect of the different shape factors $M$ and mass suction ${S}_{A}$ on the dimensionless temperature against the pseudo-similarity variable $\xi$. With an increase in ${S}_{A}$ and $M$, the solution in the upper solution branch decreases, while the solution increases in branch of second. Furthermore, solutions in the lower branch for the temperature are higher for the shape factor of lamina $M = 16.1576$ and minimum for the shape factor of the cylinder $M = 4.9$, while the reverse behavior is noticed in the upper branch. More precisely, the temperature profile (see Figure 3), the individual impact of shape factors can upsurge the TCN and improve the heat transfer rate which gives confidence that the profiles of temperature are developed higher.

Figures 4 and 5 exhibit the consequences of the parameter ${S}_{A}$ on dimensionless velocity and shear stress $f\text{'}\text{'}\left(\xi \right)$ in the entire boundary layer region of the (TiO₂-CuO/water) hybrid nanofluid for both the first and second solution branches, respectively. From Figure 4, it is observed that the first branch solution is decreasing owing to the higher value of ${S}_{A}$, while the second branch solution shows the monotonically increasing behavior. By physical observation, by the continuous rising of the impact of the mass suction, as a consequence, the hybrid nanoparticles stick with the surface of the sheet and in response, the velocity and momentum boundary layer thickness decreases. Moreover, the two distinct solution branches will not asymptotically converge at the condition far away or infinity ($\xi \to \infty$) because this condition is not the case of the considered appropriate boundary condition (3.13). The normalized velocity profiles of the hybrid nanofluid rise for the result of the first branch and decrease for the second branch due to the larger values of ${S}_{A}$ as highlighted in Figure 5. The normalized momentum boundary layer is thicker in the region approximately from $0$ to $1.6$ and then thinner in the solution region for $\xi > 1.6$. Also, the first and second solution branches hold the convergence criterion at infinity (see Figure 5) for the normalized velocity profiles.

The impact of $W{e}_{A}$ on $f\text{'}\left(\xi \right)$ and $f\text{'}\text{'}\left(\xi \right)$ distributed over the entire spatial domain for both solutions (first and second) branch versus the pseudo-similarity variable $\xi$ is portrayed in respective Figures 6 and 7. From these images, it is remarked that the normalized velocity profiles in the branch of the first solution are decelerating and escalating in the branch of the second solution as $W{e}_{A}$ increases, while the velocity profile is shrinking abruptly in both branches. Further, the thickness of the momentum boundary layer and velocity profiles are shrunken due to the rise in the value of the parameter $W{e}_{A}$ produces the resistance practiced via the hybrid nanofluid which declines its motion as seen graphically in Figure 6. Because the Weissenberg number is the ratio of relaxation time to specific process time, a decline in the specific time of the process will enhance the Weissenberg number, indicating a reduction in the velocity field. Figure 7 shows that the profile of $f\text{'}\text{'}\left(\xi \right)$ uplifts due to $W{e}_{A}$ in the first branch solution, however, the profile reduces for the branch of second results. Furthermore, the chunkiness of the normalized momentum BL is quite bigger in the outcomes of the first branch compared with the outcomes of the second branch.

Figure 8 shows the consequences of the nanoparticles volumetric fractions ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$ on the velocity profile. The outcome indicates that when the values of ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$ augment, the outcomes in the first branch show a decreasing behavior while an escalating behavior is noticed for the second branch. In general, the reduction in the profiles of velocity is observed due to the reason that when we grow up the consequence of the nanoparticles volumetric fractions produce the higher nanoparticle resistance at the wall surface of the sheet. Despite this, the velocity and momentum boundary layer shrinks owing to the higher impact of ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$.

Figures 9 and 10 elaborate the effect of ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$ on shear stress $f\text{'}\text{'}\left(\xi \right)$ distributed over the entire spatial domain and temperature profiles $\theta \left(\xi \right)$ for both distinct branch consequences, respectively. From Figure 9, it is seen that the normalized velocity profile upsurges in the first branch and declines in the second branch owing to the augmentation in ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$. In first branch outcomes, the gap in each solution curve for the corresponding varying value of ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$ is comparatively higher than the non-zoom solution of the second branch. For a better understanding of the second branch outcomes to present the gap of the solution between curves, we have shown the zooming small window as highlighted in Figure 9.

The temperature profiles for both outcome branches in response to the significant impact of ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$ are shown in Figure 10. The temperature profiles and volumetric fractions of nanoparticles hold the law of direct relation. Results reveal that both solution branches are showing an upsurge behavior owing to the higher impact of ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$. Physically, the augmentation in the thermal thickness of the boundary layer and temperature profiles is because of the improvement in the nanoparticles volumetric fractions to boost up the thermal conductivity compared to the regular fluid. Moreover, the thermal conductivity is expanding throughout the hybrid nanofluid because of the higher nanoparticles volumetric fractions, which consequences in a larger rate of heat transfer as a response to the profiles of temperature augment.

The influence of the radiation parameter ${R}_{Ad}$ on the temperature profiles of the (TiO₂-CuO/water) hybrid nanofluid for the two distinct solution branches is revealed in Figure 11. From the given figure, it is concluded that with the increasing value of ${R}_{Ad}$, the temperature profiles in both branches behave similarly like Figure 10. Generally, the heat rate of flux at the surface hastens to owe to the impact of the radiation parameter which improves the thermal conductivity as a consequence, the temperature profiles are also developed in the boundary layer region. So the thermal boundary layer thickness is improved owing to the higher impact of ${R}_{Ad}$.

Figures 12 and 13 exemplify the consequences of $W{e}_{A}$ on ${\mathit{Re}}_{x}^{1/2}{C}_{f}$ and ${\mathit{Re}}_{x}^{-1/2}N{u}_{x}$ of the (TiO₂-CuO/water) hybrid nanofluid for branch of both dissimilar solutions, respectively. In these graphs, dual solutions are obtained for the specific range of shrinking/stretching parameter, where the first solution branch exists in the whole spatial domain and the second solution survives in the range of $-\infty < {\lambda }_{A}\le -0.5$. From the graphs, it is clear that when the first and second solution branches meet at a point called CP and it is characterized by the small black solid ball but mathematically it is demarcated as ${\lambda }_{{A}_{C}}$. Therefore, the outcome of the skin friction coefficient and heat transfer exist in the range of ${\lambda }_{{A}_{C}} < {\lambda }_{A} < \infty$, no outcome is possible in the range of $-\infty < {\lambda }_{A} < {\lambda }_{{A}_{C}}$, and the unique outcome is obtained at ${\lambda }_{{A}_{C}} = {\lambda }_{A}$. Moreover, the results indicate that ${\mathit{Re}}_{x}^{1/2}{C}_{f}$ is higher in the range of $-1.5\le {\lambda }_{A} < \infty$ for the first branch solution and lower in the range of ${\lambda }_{{A}_{C}} < {\lambda }_{A} < -1.5$ owing to the higher impact of $W{e}_{A}$, while ${\mathit{Re}}_{x}^{-1/2}N{u}_{x}$ in the respective branch is decaying. In the lower branch or second solution, ${\mathit{Re}}_{x}^{1/2}{C}_{f}$ is showing a similar behavior like the first branch in range of $-\infty < {\lambda }_{A}\le -1.5$ and then starts a decomposing behavior as $W{e}_{A}$ increases while the ${\mathit{Re}}_{x}^{-1/2}N{u}_{x}$ is increasing. Also, the bifurcation values ${\lambda }_{{A}_{Ci}}$ (where $i = \mathrm{1, 2}, 3$) for each different value of the Weissenberg number. A rise in the value of the Weissenberg number causes a reduction in the value of $\left|{\lambda }_{{A}_{C}}\right|$. This tendency specifies that the higher inclusion of $W{e}_{A}$ accelerates the separation of the boundary-layer.

The variations of ${S}_{A}$ on the shear stress ${\mathit{Re}}_{x}^{1/2}{C}_{f}$ and ${\mathit{Re}}_{x}^{-1/2}N{u}_{x}$ against ${\lambda }_{A}$ are illustrated in Figures 14 and 15, respectively. There is understood a strong influence in the whole domain of the parameter ${S}_{A}$ on the shear stress and heat transfer, however, these upsurge with ${S}_{A}$ for UBS. In contrast, the skin friction coefficient is decreasing for SS owing to the rising values of ${S}_{A}$, while the heat transfer is augmented abruptly in the lower or second branch solution. In addition, it is noticed in both figures that the lower solution branch is completely terminated at ${\lambda }_{A} = -0.5$ while the first solution is existing in the whole domain. By increasing the impact of the parameter ${S}_{A}$ on the engineering gradients causes that the magnitude of the critical value $\left|{\lambda }_{{A}_{C}}\right|$ increases. In this regard, the deceleration pattern is observed in separation with the higher impact of ${S}_{A}$.

Figures 16 and 17 elucidate the influence of ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$ on the shear stress and heat transfer of the hybrid nanofluid, respectively. For growing values of ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$, the shear stress noticeably increases in FS and SB outcomes while the heat transfer increases and decreases in the first and second solution branches, respectively. Generally, the nanoparticles volumetric fractions on shear stress have an inverse relation with the velocity profiles. If we improve the impact of the hybrid nanoparticles the velocity profile is decreasing (see Figure 8) and as a result, the shear stress increases. Further, it is observed that the values of ${\varphi }_{CuO}$ and ${\varphi }_{Ti{O}_{2}}$ lead to exponential decay in the numerical critical values as mentioned in the graph windows. On the other hand, the magnitude of the bifurcation or critical value $\left|{\lambda }_{{A}_{C}}\right|$ augments. The gap in the outcomes of the first branch is emerging from the branch of the second outcomes as seen in Figures 16 and 17, respectively.

6.   Conclusions

Mathematical investigation on heat transfer and cross-flow of Williamson fluid over a porous shrinking or stretching sheet taking into account radiation and shape factors conveying hybrid nanofluid were validated in the current research. The impact of assorted pertaining parameters, like suction, volume fraction of nanoparticle, shape factor, radiation, and Williamson parameter was scrutinized. The results illustrate that double solutions are feasible through the HN which can be explained by employing certain influential parameters.

• Multiple solutions appear only when the sheet shrinks to a certain range.

• The accumulation of nanoparticle volume fraction uplifted the friction factor and improved the efficiency of the heat transfer rate of hybrid nanofluid.

• The heat transfer shrinks due to the inclusion of the Williamson parameter, while the friction factor increases.

• The temperature was enhanced due to radiation in both solutions.

• Finally, the effect of suction permits the friction factor and heat transfer to improve significantly at the porous shrinking/stretching sheet of hybrid nanofluid.

Acknowledgments

The authors are greatly acknowledged for the financial support received through the Researchers Supporting Project Number (RSP-2021/33) at King Saud University, Riyadh, Saudi Arabia.

Conflict of interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Appendix