Two-Phase Flow of Eyring–Powell Fluid with Temperature Dependent Viscosity over a Vertical Stretching Sheet

Abstract

In this work, the mixed convection flow of non-Newtonian Eyring–Powell fluid with the effects of temperature dependent viscosity (TDV) were studied together with the interaction of dust particles under the influence of Newtonian Heating (NH) boundary condition, which assume to move over a vertical stretching sheet. Alternatively, the dusty fluid model was categorized as a two-phase flow that consists of phases of fluid and dust. Through the use of similarity transformations, governing equations of fluid and dust phases are reduced into ordinary differential equations (ODE), then solved by efficient numerical Keller–box method. Numerical solution and asymptotic results for limiting cases will be presented to investigate how the flow develops at the leading edge and its end behaviour. Comparison with the published outputs in literature evidence verified the precision of the present results. Graphical diagrams presenting velocity and temperature profiles (fluid and dust) were conversed for different influential parameters. The effects of skin friction and heat transfer rate were also evaluated. The discovery indicates that the presence of the dust particles have an effect on the fluid motion, which led to a deceleration in the fluid transference. The present flow model can match to the single phase fluid cases if the fluid particle interaction parameter is ignored. The fluid velocity and temperature distributions are always higher than dust particles, besides, the opposite trend between both phases is noticed with $\beta$. Meanwhile, both phases share the similar trend in conjunction with the rest factors. Almost all of the temperature profiles are not showing a significant change, since the viscosity of fluid is high, which can be perceived in the figures. Furthermore, the present study extends some theoretical knowledge of two-phase flow.

1. Introduction

Research and studies in the area of heat and mass transport of fluids flow have discovered strategies for their development, as well as key problems. Nevertheless, it is important to choose the form of ideal fluid from the point of view of homogeneous or inhomogeneous, compressible or incompressible, Newtonian or non-Newtonian, and monophasic or polyphasic fluids, which have a significant role in determining suitable solutions for heat transfer and fluid flow enhancement. In the past few years, the new environment of accelerated technical progress has contributed to the emergence of creative approaches to analyse the suspension of fluid particles in a fluid flow that is also known as a two-phase flow model, which explains the actions of fluid dust characteristics. Industrial applications, such as petroleum transport, wastewater treatment, vehicle smoke emissions, power plant piping, and corrosive particulate matter in mining, generally involve fluid dust movement activities [1]. The movement of dust particles in a fluid leads to a two-phase cycle. Numerous pieces of research on the dynamics of solid particles in the fluid have been conducted as a result of recent advancements in the field of two-phase flow. The fluid and solid phases of this solid–liquid system are independently formulated using different continuum equations. This phenomenon includes micro-propulsion, aerosol filtration, powder transport, oil industry, and the flow of corpuscle in plasma (a liquid with suspended solids). Overall, it is very useful for modelling flow with a binary mixture of non-Newtonian fluid and solid particles linked to certain conditions. It can thus be suggested that this two-phase model could benefit in studying the dusty Eyring–Powell fluid that exhibits the binary characteristics of the Eyring–Powell fluid and spherical dust particles, such as undertaken here. Ref. [2] studied the process of radiative heat transfer in the flow of dusty liquid under the power generation aspects. The boundary layer flow of a dusty fluid with electrically conducting criteria in a porous medium has been studied by [3]. When the interaction of these phases is significant, the temperature of the fluid is always higher than that of the dust. In accordance with these applications, the literature includes a variety of works in corresponding flow for various contexts, such as multiple geometries, boundary conditions and fluid-based forms. The non-Newtonian Casson model with dust particles have been developed by [4,5]. In addition, [6,7] utilized the treating fluid–particle interaction with buoyancy forces on Jeffrey fluid with Newtonian heating indicates that the presence of the dust particles has an effect on the fluid motion, which led to decelerate the fluid transference. The natural convection flow caused by non-Newtonian fluid with dust nanoparticles has been addressed in [8,9]. Furthermore, [10] theoretically analysed along a vertical stretching sheet for magnetohydrodynamic (MHD) mixed convection of non-Newtonian tangent hyperbolic nanofluid flow with suspended dust particles. In [11], a detailed study has been done on a two-phase model implemented in the presence of hybrid nanoparticles on the dusty liquid flow through a stretching cylinder by employing the modified Fourier heat flux law. The study was conducted by considering the effect of viscous dissipation and non-linear thermal radiation, which demonstrated a two-phase dusty liquid movement across a permeable surface. Other contributions of flow models on dusty non-Newtonian fluid have drawn substantial interest among researchers under different conditions [12,13,14,15].

Consideration of the study on boundary layer phenomenon focussing on non-Newtonian heat-transported substances is essential for a deeper comprehension of engineering and industrial–technology issues, and the movement of these materials occurs extensively in various industrial processes, such as guided missiles, rain erosion, fluidisation, atmospheric failure, lunar ash fall, paint and aerosol spraying, as well as the cooling of nuclear reactors. Although the existence of such substances is greatly complicated and troublesome, a variety of constitutive models have been developed and studied to research the correct flow behaviour. Eyring–Powell fluid model is one of the subcategories of the non-Newtonian fluid model. It has a clear characteristic under other non-Newtonian models, conveniently derived from the kinetic theory of gases rather than empirical relations and comes baked from Newtonian behaviour for low and high shear rates. The rheological paradigm is known for its robustness and versatility in physical action. Ref. [12] addressed movement attributable to pulsatile pressure gradient of dusty non-Newtonian fluid with heat transfer in a channel. A preliminary analysis of the magnetohydrodynamic movement of the Eyring–Powell liquid under the suspension of nanoparticles and dust has been done by [13] in which this model has shown that the intensity of heat transfer in the aluminium oxide nanofluid was higher than that in the ferro oxide nanofluid with the current viscous variance parameter. In [14], by considering variable thermal conductivity and thermal radiation, analytical solutions of unstable flow Eyring–Powell and Carreau non-Newtonian fluids in the suspension of dust and nickel nanoparticles, a higher heat transfer rate was recorded in the nickel + Eyring–Powell mixture compared to the nickel + Carreau case. Ref. [15] attempted to examine the effects of heat and mass transfer in the presence of nonlinear convection and thermal radiation of MHD rheological Eyring–Powell fluid with dust and graphene nanoparticles in a mixture of ethylene glycol. In order to study the relationship between the fluid and dust phased, this analysis was conducted to further investigate the flow behaviour of dusty Powell–Eyring fluid in the vertical stretching sheet associated with temperature-dependent viscosity combined with NH as a thermal boundary conditions.

Over recent decades, the thermal boundary conditions have a great influence on the heat transfer in the laminar boundary layer flow problem. However, there is a situation in which the heat transfer rate is proportional to the local difference in temperature with ambient conditions or usually termed conjugate boundary conditions, which are driven by NH, should be considered as well. The case of NH has been mentioned by [16] during a study on the boundary layer flow over an upright plate. Ref. [17] made further headway during an investigation on free convection flow across horizontal surface. Ref. [18] proposed a series of solutions and numerical Eyring fluid flow with NH. Ref. [19] provided detailed nonlinear convective magneto nanofluid Eyring fluid with the effects of NH, while Ref. [20] included the effects of the thermal radiation Eyring fluid subjected to NH boundary condition. Meanwhile, Refs. [21,22] highlighted the effects of MHD and NH on Powell–Eyring fluid over a stretching cylinder and inclined permeable surface, respectively. Ref. [23] reviewed the problem of temperature-dependent viscosity on mixed convection flow of Eyring–Powell fluid studied together with NH. A mathematical model of forced convective flow on non-Newtonian Eyring–Powell fluid under temperature-dependent viscosity circumstance is formulated by [24]. The effects of non-Newtonian magnetohydrodynamic nanofluid over a stretched plate with NH effect have been investigated by [25]. Other studies associated with particular impact for fluid–solid flow considering different fluid models were established and reported in [26,27,28,29].

Motivated by the impactful research as scrutinised above, this present study is dedicated to examining the two-phase boundary layer flow of Eyring–Powell fluid together with the temperature-dependent viscosity from a vertical stretching sheet, the temperature of which is higher than that of the ambient fluid. The simulation of mixed convection influence with NH was also implanted in this investigation. The mathematical formulations of the problems are constructed as mentioned in the study scope, which involves the derivation of governing equations for the proposed problem. The similarity transformation is used to transform the non-linear governing equations into ordinary differential equations (ODE). Then, the numerical solutions of the transformed equations are solved using the Keller-box method.

The step size of time and space can be arbitrary, since this method is implicit with second order accuracy, which makes it suitable to solve the parabolic partial differential equations efficiently [30]. However, the computation could be time consuming if the small step size of time and space is inserted. The algorithm of the Keller-box method is computed in Matlab software to generate the results and the figures for various non-dimensional parameters on the velocity and temperature profile. The comparison with the previous published result were tabulated to verify the present results by fixing several parameters. The output from the investigation are useful for the scientist and experimentalist in studying the behaviour of fluids, which have interactions with dust particles.

2. Mathematical Formulation

Flow suspended with particles affected by TDV over a vertical stretching sheet was introduced for the steady incompressible mixed convection of non-Newtonian Eyring fluid under NH condition. The term, $T_{\infty }$, is related to temperature of ambient fluid. The $x$-axis was oriented to the vertical plane, and the $y$-axis to the plane was perpendicular. The sheet is stretched with the velocity, $u_w(x)=bx$ where, $b>0$ is stretching rate. The flow was created by the stretching of the sheet due to the simultaneous application of two equal and opposite forces along the $x$-axis, holding the origin fixed and finding the flow to be limited to the area, $b>0$. The configuration of a physical model is displayed in Figure 1.

The solid particles are presumed to be spherical in shape and uniform size, where their density remains constant and the inter-particle collision may be neglected, since they are considered to be diluted throughout the flow. In terms of heat transfer from one atom to another, all of these were ignored: the volume fraction of dust particles, coagulation, phase transition and deposition. The fluid and sediment particle motions were linked only through drag-and-heat flow across them, the stokes linear drag theory was used for modelling the drag force. Based on the preceding assumptions, the basic two-dimensional boundary layer equations involving continuity, momentum and energy for both the fluid and particle phases with usual ratings can be written as [31,32] shown below:

For fluid phase:

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

(1)

$$u\frac{\partial u}{\partial x}+v\frac{\partial u}{\partial y}=\frac{1}{\rho }\frac{\partial }{\partial y}\left(\mu \frac{\partial u}{\partial y}\right)+\frac{1}{\rho \tilde{\beta }c}\left(\frac{{\partial }^{2}u}{\partial y^2}\right)-\frac{1}{2\rho \tilde{\beta }{c}^3}{\left(\frac{\partial u}{\partial y}\right)}^2\frac{{\partial }^{2}u}{\partial y^2}+\beta g\left(T-T_{\infty }\right)+\frac{{\rho }_p}{\rho \tau }\left(u_p-u\right),$$

(2)

$$\rho c_p\left(u\frac{\partial T}{\partial x}+v\frac{\partial T}{\partial y}\right)=k\left(\frac{{\partial }^{2}T}{\partial y^2}\right)+\frac{{\rho }_p{c}_s}{{\gamma }_T}\left(T_p-T\right).$$

(3)

For dust phase:

$$\frac{\partial u_p}{\partial x}+\frac{\partial v_p}{\partial y}=0,$$

(4)

$${\rho }_p\left(u_p\frac{\partial u_p}{\partial x}+v_p\frac{\partial u_p}{\partial y}\right)=\frac{{\rho }_p}{{\tau }_v}(u-u_p),$$

(5)

$${\rho }_p{c}_s\left(u_p\frac{\partial T_p}{\partial x}+v_p\frac{\partial T_p}{\partial y}\right)=-\frac{{\rho }_p{c}_s}{{\gamma }_T}(T_p-T).$$

(6)

Here, $\left(u,v\right), T, \rho , c_p$ and $\mu$ represents the components of velocity in $\left(x, y\right)$ directions, temperature, density, specific heat at constant pressure and viscosity coefficient, respectively. Meanwhile, $\left(u_p,v_p\right),T_p, {\rho }_p, c_s, {\tau }_v$ and ${\gamma }_T$ denote the velocity components in $\left(x,y\right)$ directions, temperature, density, specific heat, velocity and thermal relaxation time for dust phase, respectively. The corresponding fluid and particle phase boundary conditions were given as.

$$\begin{array}{l}u=u_w(x)=ax, v=0, \frac{\partial T}{\partial y}=-h_{s}T \mathrm{at} y=0\\ u\to 0, u_p\to 0, v_p\to v, T\to T_{\infty }, T_p\to T_{\infty } \mathrm{as} y\to \infty \end{array}$$

(7)

In (7), the parameter was corresponded to the velocity of the stretching surface with $a$ being a positive constant of $u_w(x)$, thermal conductivity $k$, heat transfer coefficient $h_s$ and ambient temperature $T_{\infty }$. To obtain the set of similarity equations in the form of ordinary differential equations, the similarity transformations (8) were adopted and applied to the governing Equations (1)–(6).

$$\begin{array}{l}u=ax{f}^{\prime }(\eta ), v=-{\left(a\upsilon \right)}^{1/2}f(\eta ), \eta ={\left(\frac{a}{\upsilon }\right)}^{1/2}y,\\ \theta (\eta )=\frac{T-T_{\infty }}{T_{\infty }}, u_p=ax{F}^{\prime }(\eta ), v_p=-{\left(a\upsilon \right)}^{1/2}F(\eta ), {\theta }_p(\eta )=\frac{T_p-T_{\infty }}{T_{\infty }},\end{array}$$

(8)

The Reynolds exponential viscosity model was used to predict temperature-dependent variability in viscosity that gave a detailed approach as:

$$\mu \left(\theta \right)={\mu }_0{e}^{-\left({\beta }_1\theta \right)}={\mu }_0\left[1-\left({\beta }_1\theta \right)+O\left({\beta }_1{}^2\right)\right],$$

(9)

The Equations (1)–(6) are changeable from PDEs to ODEs, which can be represented as:

$$\begin{array}{l}\left(1+M\right)f^{‴}(\eta )-{\left(f^{\prime }\left(\eta \right)\right)}^2+f(\eta )f^{″}\left(\eta \right)+\beta N\left(F^{\prime }(\eta )-f^{\prime }(\eta )\right)\\ -BM{\left(f^{″}\left(\eta \right)\right)}^2{f}^{‴}\left(\eta \right)-\alpha f^{″}\left(\eta \right){\theta }^{\prime }\left(\eta \right)-\alpha \theta \left(\eta \right)f^{‴}(\eta )+\lambda \theta =0,\end{array}$$

(10)

$${\theta }^{″}\left(\eta \right)+\mathrm{Pr}f(\eta ){\theta }^{\prime }(\eta )+\frac{2}{3}\beta N\left({\theta }_p(\eta )-\theta (\eta )\right)=0,$$

(11)

$${\left(F^{\prime }(\eta )\right)}^2-F(\eta )F^{″}(\eta )+\beta \left(F^{\prime }(\eta )-f^{\prime }(\eta )\right)=0,$$

(12)

$${\theta }_p{}^{\prime }(\eta )F(\eta )+\frac{2}{3}\frac{\beta }{\mathrm{Pr}\gamma }\left(\theta (\eta )-{\theta }_p(\eta )\right)=0$$

(13)

Subjected to boundary conditions:

$$\begin{array}{l}f(0)=0,f^{\prime }(0)=1, {\theta }^{\prime }(0)=-{\gamma }_1\left(1+\theta (0)\right) \mathrm{at} \eta =0\\ f^{\prime }(\eta )\to 0, F^{\prime }(\eta )\to 0, F(\eta )\to f(\eta ),\\ \theta (\eta )\to 0, {\theta }_p(\eta )\to 0 \mathrm{as} \eta \to \infty \end{array}$$

(14)

In Equations (10)–(14), a notation prime $({}^{\prime })$ corresponds to the differentiation with respect to $\eta .$ Additionally, the dimensionless numbers and parameters are as follows where $M=\frac{1}{{\mu }_0\tilde{\beta }c}$ and $B=\frac{a^3{x}^2}{2c^2{\nu }_f}$ are the fluid parameters, $\mathrm{Pr}={\nu }_f/\alpha$ represents Prandtl number, $\alpha =k/\rho c_p$ represents viscosity parameter, ${\gamma }_1$ conjugate parameter of heat transfer, $\lambda =g_c{\beta }_T(T_f-T_{\infty })/a^{2}x\lambda$ is mixed convection parameter, $\gamma =c_s/c_p$ specific heat ratio of mixture parameter, $N={\rho }_p/\rho$ parameter of mass concentration of particle phase, fluid–particle interaction parameter $\beta =1/a{\tau }_v$ and Reynolds number ${\mathrm{Re}}_x=(a{x}^2/\nu )$. A limiting case arising in this problem was without the presence of dust particles effect where the buoyancy force is negligible and can be obtained using the following expression [33]:

$$f(\eta )=S+A\eta +\left(1-A\right)\left(1-\exp \left(-\eta \right)\right), \theta (\eta )=\frac{{\gamma }_1}{1-{\gamma }_1}\exp \left(-\eta \right).$$

(15)

It is important to mention here that the comparison between the present results with the exact solution is necessary to claim the accuracy of the current model and its output. The primary physical quantity of importance is the dimensionless coefficient of skin friction and the local Nusselt number, which has been described by (16), where the shear stress and surface heat are compatible with those referred in [34].

$$C_{fx}=\frac{{\tau }_w}{\rho U_w^2(x)}, N{u}_x=\frac{x{q}_w}{k\left(T_w-T_{\infty }\right)}$$

(16)

where

$${\tau }_w=\left({\mu }_0+\frac{1}{\tilde{\beta }{C}^{*}}\right)\frac{\partial u}{\partial y}-\frac{1}{6\tilde{\beta }}{\left(\frac{1}{C^{*}}\frac{\partial u}{\partial y}\right)}^3, \mathrm{and} q_w=-k{\left(\frac{\partial T}{\partial y}\right)}_{y=0}$$

(17)

The shear stress and surface heat transfer is calculated using the following definition:

$$C_f{\mathrm{Re}}_x^{1/2}=\left(\left(1-\alpha \theta \left(0\right)\right)+M\right)f^{″}(0)-\frac{B}{3}M{f^{″}}^3\left(0\right), N{u}_x{\mathrm{Re}}_x{}^{-1/2}={\gamma }_1\left(\frac{1}{\theta (0)}+1\right)$$

(18)

3. Results and Discussion

Equations (10)–(13) were solved numerically, along with boundary conditions (14) using the Keller-box approach as computed in the Matlab program. The Keller-box method comprise of four steps which are:

Step 1: The nonlinear partial differential equation are first transformed to first order system;

Step 2: The first order system is then approximated using central difference;

Step 3: The Newton’s method is applied to linearize the system;

Step 4: The linearized system is solved by block elimination technique;

To initially integrate the procedure, the nonlinear structure of ordinary differential equations was converted into a structure of linear first-order equations. Our bulk calculations were viewed with ${\eta }_{\infty }=8$ and identified as appropriate for all values of the parameters, which were considered asymptotically to achieve the far-field boundary conditions as seen 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, Figure 15, Figure 16 and Figure 17. The interaction force among two phases is significant, where usually, both governing equations are coupled through the term of total fluid–particle interaction force per unit volume that is clearly different from single phase flow. Numerical computation are conducted for fluid parameters, namely $M$ and $B$, Prandtl number $\mathrm{Pr}$, viscosity parameter $\alpha$, mixed convection parameter $\lambda$, ${\gamma }_1$ conjugate parameter of heat transfer, a parameter of mass concentration of particle phase $N$ and fluid–particle interaction parameter $\beta$.

The local Nusselt number is one of the important characteristics in the heat transfer field, which indicates the ratio of convective heat transfer to conductive heat transfer. Hence, in order to check the accuracy of the numerical method used, the comparison in Nusselt number $N{u}_x{\mathrm{Re}}_x{}^{-1/2}$ for a fixed value of $\mathrm{Pr}$ with the established results in [35,36,37] and it is revealed to be in strong agreement as displayed in

Table 1. Comparative study on value $-{\theta }^{\prime }\left(0\right)$.

$\mathbf{Pr}$	[35]	[36]	[37]	Present
1	1.3333	1.3333	1.3333	1.3329
3	2.50970	2.50972	2.50972	2.50969
10	4.79690	4.79686	4.79687	4.79689
100	15.7120	15.7118	15.7120	15.7098

. In addition, a direct comparative study was carried out with the exact Equation (15), as well as the existing study reported with the available published result by [38,39,40,41] shown in

Table 2. Comparative study on $f^{″}\left(0\right)$.

Existing Literature	Model of Problem	Boundary Condition	Limiting Cases	$\mathbf{Value}\mathbf{of}{\mathit{f}}^{″}\left(0\right)$
Exact Solution (15)	$f^{″}=-e^{-\eta }$	-	-	−1.0000
[38]	$f^{‴}-{f^{\prime }}^2+f{f}^{″}-A\left(f^{\prime }+\frac{1}{2}\eta f^{″}\right)+\lambda \theta =0$	$\begin{array}{l}f\left(0\right)=0\\ f^{\prime }\left(0\right)=1\\ f^{\prime }(\infty )=0\end{array}$	$A=\lambda =0$	−1.0000
[39]	$f^{‴}-{f^{\prime }}^2+f{f}^{″}-k_1\left(2f^{\prime }{f}^{‴}-{f^{″}}^2-f{f}^{‴}\right)=0$	$\begin{array}{l}f\left(0\right)=0\\ f^{\prime }\left(0\right)=1\\ f^{\prime }(\infty )=0\end{array}$	$k_1=0$	−1.0000
[40]	$\begin{array}{l}\left(1+M\right)\left(1+2\eta \gamma \right)f^{‴}-\alpha M{\left(1+2\eta \gamma \right)}^2{f^{″}}^2{f}^{‴}\\ +2\gamma \left(1+M\right)f^{″}+f{f}^{″}-\frac{4}{3}\alpha M\left(\gamma +2\eta {\gamma }^2\right){f^{″}}^3\\ -{f^{\prime }}^2+\lambda \theta \sin \phi =0\end{array}$	$\begin{array}{l}f\left(0\right)=0\\ f^{\prime }\left(0\right)=1\\ f^{\prime }(\infty )=0\end{array}$	$\begin{array}{l}\alpha =\gamma =0\\ \lambda =M=0\end{array}$	−1.0000
[41]	$\left(1+M\right)f^{‴}-{f^{\prime }}^2+f{f}^{″}-MB{f^{″}}^2{f}^{‴}-H{f}^{\prime }=0$	$\begin{array}{l}f\left(0\right)=0\\ f^{\prime }\left(0\right)=1\\ f^{\prime }(\infty )=0\end{array}$	$\begin{array}{l}B=H=0,\\ M=0.0001\end{array}$	−1.0000
Present study	$\begin{array}{l}\left(1+M\right)f^{‴}-{f^{\prime }}^2+f{f}^{″}+\beta N\left(F^{\prime }-f^{\prime }\right)\\ -BM{f^{″}}^2{f}^{‴}-\alpha f^{″}{\theta }^{\prime }-\alpha \theta f^{‴}+\lambda \theta =0,\end{array}$	$\begin{array}{l}f\left(0\right)=0\\ f^{\prime }\left(0\right)=1\\ f^{\prime }(\infty )=0\end{array}$	$\begin{array}{l}B=M=0\\ N=\beta =0\\ \alpha =\lambda =0\end{array}$	−1.0015

. From Table 1 and Table 2, an excellent agreement is achieved, which indicates that the current model and its findings are acceptable. It is worth declaring here, even in the limiting cases, the present model does not exactly give the same solution, but the difference is very small. It is logical since the present model is more complex with multiple parameters.

Table 3. Numerical results of $C_f{\mathrm{Re}}_x{}^{1/2}$ and $N{u}_x{\mathrm{Re}}_x{}^{-1/2}$ for various values of $\mathrm{Pr}, \alpha , M, B, \lambda , {\gamma }_1, \beta$ and $N$.

$\mathbf{P}\mathbf{r}$	$\mathit{\alpha }$	$\mathit{M}$	$\mathit{B}$	$\mathit{\lambda }$	${\mathit{\gamma }}_1$	$\mathit{\beta }$	$\mathit{N}$	${\mathit{C}}_{\mathit{f}}{\mathbf{Re}}_{\mathit{x}}{}^{1/2}$	$\mathit{N}{\mathit{u}}_{\mathit{x}}{\mathbf{Re}}_{\mathit{x}}{}^{-1/2}$
7	0.1	0.6	0.6	0.1	0.5	0.5	0.5	−1.076579	0.105283
9								−1.078272	0.104603
12								−1.079770	0.103939
10	0.2							−1.077788	0.104348
	0.4							−1.075612	0.104349
	0.7							−1.072307	0.104350
10	0.1	0.5						−1.048326	0.104357
		0.9						−1.173040	0.104324
		1.5						−1.351057	0.104292
10	0.1	0.6	0.1					−1.133436	0.104328
			0.5					−1.088912	0.104324
			0.9					−1.048210	0.104322
10	0.1	0.1	0.6	0.3				−1.754138	0.342542
				0.5				−1.741401	0.630439
				0.9				−1.706099	1.433335
10	0.1	0.6	0.6	0.1	0.1			−1.085300	0.104326
					0.6			−1.081074	0.104324
					1.2			−1.076098	0.104318
10	0.1	0.6	0.6	0.1	0.5	0.1		−1.180434	0.104421
						0.4		−1.114385	0.104359
						0.9		−1.044242	0.104282
10	0.1	0.6	0.6	0.1	0.5	0.5	0.1	−1.182364	0.104415
							0.4	−1.125083	0.104364
							0.9	−1.019745	0.104280

demonstrates the variance of the skin friction coefficient and the Nusselt numbers for various parameters of the present analysis.

Figure 2, Figure 3, Figure 4 and Figure 5 was plotted to understand the velocity and temperature distribution of the fluid and particle phase under variance of $\mathrm{Pr}$ and $\alpha$. It was revealed that the velocity were decreased for both fluid and particle phases as $\mathrm{Pr}$ and $\alpha$ increases. The similar trend was noticed in temperature distribution for both phases in increasing $\mathrm{Pr}$ but contrary in the growing of $\alpha$. At far from the surface, it is remarked the profile asymptotically reached the boundary conditions, and, therefore, the authors are confident on the correctness of present results.

Figure 6, Figure 7, Figure 8 and Figure 9 display the distribution on velocity and temperature for multiple values of $M$ and $B$. It was found that with higher elasticity parameter (presence Eyring fluid), the magnitude of velocity for both fluid and particle were enhanced. The change in the velocity contributed to boosting the heat of the fluid (for increasing $M$) but against the heating development (for increasing $B$).

Figure 10 and Figure 11 indicate that the amount of ${\gamma }_1$. Physically, the heat transfer rate with the ${\gamma }_1$ was reduced; this reduced the temperature and resulting thickness of the boundary layer. During the event with large values of the heat transfer equation, which implies a strong heat transfer rate as it slowly declines, the fluid became a small heat transfer rate. Through mixed convection and the effects of the fluid–particle interaction parameters $\beta$ and $\lambda$, the velocity of all the phases and their associated boundary layer thicknesses increased. The response to the temperature profile was quite the opposite of the speed field with $\lambda$ and $\beta$, as seen in Figure 12, Figure 13, Figure 14 and Figure 15, further finding that the flow properties of the dusty fluid can be greatly regulated by changing the influence of the parameters of fluid–particle interaction $\beta$.

Figure 16 and Figure 17 were plotted to evaluate all fluid and particle phase velocity and temperature components for the variance of $N$. Within the boundary layer, temperature profiles reduced with improvement within $N$. On the other hand, the velocity profile was improved by increasing the parameter $N$. This was because the fluid tends to raise the intensity of drag between the phases with the mass content of dust particles rising. The fluid movement was thus slowed down, resulting in reduced surface-phase energy, since the surface layer was pulled together with the liquid. By continuing to increase the mass content of the dust particles, more fluid-phase energy was converted into a larger number of particles, but less energy from the fluid phase was supplied to the individual particles. Therefore, it can be inferred that varying $N$ will greatly affect the flow characteristics. Furthermore, the boundary momentum layer for ordinary Eyring fluid was observed thinner than that of the dusty Eyring fluid.

Figure 2. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=\lambda =0.5$, ${\gamma }_1=\alpha =0.1$ and $\beta =N=0.6$ for various values of $\mathrm{Pr}$.

Figure 2. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=\lambda =0.5$, ${\gamma }_1=\alpha =0.1$ and $\beta =N=0.6$ for various values of $\mathrm{Pr}$.

Figure 3. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=\lambda =0.5$, ${\gamma }_1=\alpha =0.1$ and $\beta =N=0.6$ for various values of $\mathrm{Pr}$.

Figure 3. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=\lambda =0.5$, ${\gamma }_1=\alpha =0.1$ and $\beta =N=0.6$ for various values of $\mathrm{Pr}$.

Figure 4. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=\lambda =0.5,{\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $\alpha$.

Figure 4. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=\lambda =0.5,{\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $\alpha$.

Figure 5. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=\lambda =0.5,{\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $\alpha$.

Figure 5. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=\lambda =0.5,{\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $\alpha$.

Figure 6. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $B=\lambda =0.5,\alpha ={\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $M$.

Figure 6. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $B=\lambda =0.5,\alpha ={\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $M$.

Figure 7. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $B=\lambda =0.5,\alpha ={\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $M$.

Figure 7. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $B=\lambda =0.5,\alpha ={\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $M$.

Figure 8. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=\lambda =0.5,\alpha ={\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $B$.

Figure 8. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=\lambda =0.5,\alpha ={\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $B$.

Figure 9. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=\lambda =0.5,\alpha ={\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $B$.

Figure 9. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=\lambda =0.5,\alpha ={\gamma }_1=0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $B$.

Figure 10. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=\lambda =0.5,\alpha =0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of ${\gamma }_1$.

Figure 10. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=\lambda =0.5,\alpha =0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of ${\gamma }_1$.

Figure 11. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=\lambda =0.5,\alpha =0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of ${\gamma }_1$.

Figure 11. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=\lambda =0.5,\alpha =0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of ${\gamma }_1$.

Figure 12. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=0.5,{\gamma }_1=\alpha =0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $\lambda$.

Figure 12. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=0.5,{\gamma }_1=\alpha =0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $\lambda$.

Figure 13. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=0.5,{\gamma }_1=\alpha =0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $\lambda$.

Figure 13. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=0.5,{\gamma }_1=\alpha =0.1$, $\beta =N=0.6$ and $\mathrm{Pr}=10$ for various values of $\lambda$.

Figure 14. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=\lambda =0.5,N=0.6$, ${\gamma }_1=\alpha =0.1$ and $\mathrm{Pr}=10$ for various values of $\beta$.

Figure 14. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=\lambda =0.5,N=0.6$, ${\gamma }_1=\alpha =0.1$ and $\mathrm{Pr}=10$ for various values of $\beta$.

Figure 15. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=\lambda =0.5,N=0.6$, ${\gamma }_1=\alpha =0.1$ and $\mathrm{Pr}=10$ for various values of $\beta$.

Figure 15. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=\lambda =0.5,N=0.6$, ${\gamma }_1=\alpha =0.1$ and $\mathrm{Pr}=10$ for various values of $\beta$.

Figure 16. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=\lambda =0.5,\beta =0.6$, ${\gamma }_1=\alpha =0.1$ and $\mathrm{Pr}=10$ for various values of $N$.

Figure 16. $f^{\prime }\left(\eta \right)$ and $F^{\prime }\left(\eta \right)$ at $M=B=\lambda =0.5,\beta =0.6$, ${\gamma }_1=\alpha =0.1$ and $\mathrm{Pr}=10$ for various values of $N$.

Figure 17. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=\lambda =0.5,\beta =0.6$, ${\gamma }_1=\alpha =0.1$ and $\mathrm{Pr}=10$ for various values of $N$.

Figure 17. $\theta \left(\eta \right)$ and ${\theta }_p\left(\eta \right)$ at $M=B=\lambda =0.5,\beta =0.6$, ${\gamma }_1=\alpha =0.1$ and $\mathrm{Pr}=10$ for various values of $N$.

4. Conclusions and Future Work

The results of temperature-dependent viscosity on the mixed convection flow of non-Newtonian Eyring–Powell fluid due to a vertical stretching surface interacting with dust particles under the control of NH boundary conditions have been numerically explained. Highly non-linear governing PDEs were first converted to ODEs, utilising similarity transformations and then numerically resolved using the Keller-box method. MATLAB program has been applied to conduct computational. Velocity and temperature profile for all fluid and particle phases were illustrated graphically. The value of sheer rate and heat transfer coefficient were calculated and have tabulated in table for various pertinent parameters. In addition, the current output was validated by comparative study with previously reported results and perceived a good agreement between them. In conjunction, a few important findings from the existing research are indicated as:

• For certain applications, the fluid’s flow and heat transfer can be regulated by embedding the particles of fine dust.
• The process velocity (temperature) of fluid and dust have the opposite effects for buoyancy force parameter variability.
• In the mixed convection regime, the local shear stress increases and the local rate of heat transfer decreases as the value of buoyancy parameter increases for all values of the Prandtl number and the viscosity variation parameter.
• The fluid–particle interaction parameter variability is favourable for the thickness of the dust boundary layer. However, for the thickness of boundary layer of momentum, it is unfavourable.
• The velocity distribution was suppressed with the Prandtl number compared to temperature distribution.
• The velocity profiles increase and the viscosity of the fluid decrease near the surface of the plate owing to increase in the value of the viscosity variation parameter. The temperature profiles of both phases are enhanced by rising $\alpha$.
• The quantity of skin friction decreases with greater values of fluid parameters, mixed convection, conjugate parameter of heat transfer, mass concentration of particle phase and fluid–particle interaction parameter.
• Increase in the value of the viscosity variation parameter leads to increase in the local shear stress and to decrease in the local rate of heat transfer. Its effect on the increase of the rate of heat transfer is less than that of the local shear stress.

Ultimately, it is worth concluding that the existence of dust particles has a substantial effect on the flow behaviour of Eyring fluid in the presence of TDV.