Numerical Investigation of Mixed Convective Williamson Fluid Flow Over an Exponentially Stretching Permeable Curved Surface

: The present investigation aims to examine the heat ﬂux mechanism in the hagnetohydro-dynamic (MHD) mixed convective ﬂow of Williamson-type ﬂuid across an exponential stretching porous curved surface. The signiﬁcant role of thermal conductivity (variable), non-linear thermal radiation, unequal source-sink, and Joules heating is considered. The governing problems are obtained using the Navier–Stokes theory, and the appropriate similarity transformation is applied to write the partial differential equations in the form of single-variable differential equations. The solutions are obtained by using a MATLAB-based built-in bvp4c package. The vital aspect of this analysis is to observe the effects of the curvature parameter, magnetic number, suction/injection parameter, permeability parameter, Prandtl factor, Eckert factor, non-linear radiation parameter, buoyancy parameter, temperature ratio parameter, Williamson ﬂuid parameter, and thermal conductivity (variable) parameter on the velocity ﬁeld, thermal distribution, and pressure proﬁle which are discussed in detail using a graphical approach. The correlation with the literature reveals a satisfactory improvement in the existing results on permeability factors in Williamson ﬂuids.


Introduction
In modern research, non-Newtonian fluids have great significance. The rheological attributes of non-Newtonian type fluids cannot be well illustrated by the famous Navier-Stokes equation only. Therefore, scientists proposed several models to portray the characteristics of non-Newtonian fluids. In 1929, Williamson proposed a model for non-Newtonian fluids that illustrates the rheological attributes of non-Newtonian type fluids. In recent years, many scientists worked on the Williamson model to study the behavior of non-Newtonian fluids. Kebede et al. [1] recently reported the heat and mass transport in Williamson nanofluids. Kumar et al. [2] investigated the MHD fluid flow bounded by a curved sheet for the Williamson fluid model. They found increasing outcome of rising curvature parameter on velocity and decrement in it is seen for large values of the Williamson and magnetic parameters. In a study, Nadeem et al. [3] analyzed the flow bounded by a stretching sheet for Williamson type fluid model and discussed the results graphically. Nadeem et al. [4] investigated the flow of Williamson-type nanofluid determined by stretching sheet by presenting a new form of governing equations. Salahuddin et al. [5] analyzed the Induced MHD impact on exponentially varying viscosity of Williamson fluid flow with variable conductivity and diffusivity. Ibrahim [6] used the Williamson model to tackle the MHD flow over a stretchable cylinder using the porous medium. Iqbal et al. [7] recently analyzed the Williamson fluid flow in a cylinder by considering an exponentially stretching sheet. They presented the effects on velocity and skin friction graphically by adopting the shooting technique. Panezai et al. [8] numerically investigated the heat transfer characteristics of Williamson fluid over a porous Wedge. Salahuddin et al. [9] studied the boundary layer phenomena of Williamson-type fluid with slip conditions over a stretching cylinder.
Heat transfer flows with the porous surface are of considerable interest due to their advanced engineering applications and their frequent use in industrial technology and power generation systems. Food processing, cooling towers, and the distribution of temperature are the leading examples. Kairi et al. [10] reported the effect of melting on mixed convection heat and mass transfer in a non-Newtonian-fluid-saturated, non-Darcy porous medium. Bestman [11] explored the heat and mass flux accompanying the free convection in the boundary layer flow using the porous medium. Lakshmi et al. [12] reported mixed convection type stagnation point fluid flow analysis using a stretching type porous surface. Ambreen et al. [13] reported an analytic investigation on boundary layer phenomena and heat flux, specifically using a rotating disk with porous attributes. In a study, Rasool et al. [14] reported the entropy (disorder) consequences and the high impact of a second-order chemical reaction in a nanofluid based on the Darcy model known for porous medium using a non-linear stretching surface. Shafiq et al. [15] presented the importance of convective boundary and thermal slip in a three-dimensional and rotating disk (porous medium). Reddy et al. [16] presented numerical results for the boundary layer phenomena of a nanofluid using a porous medium on an exponentially stretching sheet.
Magnetohydrodynamics (MHD) is the study of fluid behavior in the presence of electric and magnetic fields. MHD plays a significant role in biomedical and industrial sciences [17]. Some applications of MHD are magnetic drug targeting, cancer tumor treatment, magnetic devices for cell separation, magnetic endoscopy, and adjusting blood flow during surgery. Much advancement is observed in research related to MHD flows over a stretching surface in the literature. Jabeen et al. [18] considered a porous surface to study the effects of MHD boundary layer flow over a stretching surface and observed that Maxwell fluid has a high thermal conductivity rate compared to other non-Newtonian fluid models. Rasool et al. [19] analyzed the high impact of Darcy medium (porous medium) on MHD nanofluid bounded by a stretching (non-linear) surface. In another study, the authors Rasool et al. [20] use the famous SPECTRAL method to figure out the impact of EMHD by a vertical Riga plate in a second-grade type nanofluid. Akbar et al. [21] numerically investigated the MHD flow of second-grade fluids in a porous medium with prescribed vorticity. Aurangzaib et al. [22] reported the unsteady MHD type mixed convective flow based on a micropolar porous surface.
In recent research, the boundary layer flow over stretching curved/flat surfaces has attracted much attention for specific reasons, such as their extensive industry and engineering sector applications. Many fundamental processes involve condensation, drawing of wires, metal extrusion, rolling, fiber spinning, and polymer sheet extrusion, showing the essential applications of stretching surfaces in industries. Sakiadis [23] derived the fundamental differential equations for the boundary layer type theory of continuous solid surfaces for both the laminar and turbulent flow. Crane [24] proposed an exact solution for surface friction and thermal conductivity for the flow past a stretching surface. Gupta et al. [25] investigated mass, heat, and momentum over a stretching surface in the boundary layer subject to blowing and suction. Okechi et al. [26] investigated the fluid flow phenomena along exponential and stretching curved surface. Kumar et al. [27] analyzed the MHD fluid flow of micropolar fluid passed over an exponential and stretching carved surface. Nayak [28] investigated MHD convective flow to study the impact of thermal diffusion via an exponentially stretching sheet and proposed an iterative solution by using the R.K. method and shooting technique. Jalil and Asghar [29] presented a pioneer work by discussing exponential stretching using Lie group analysis for the first time. They used the perturbation method to extend their investigation for shear-thinning fluids. Sanni et al. [30] blended the stretching surface and curved surface concepts to study viscous fluid flow. Recently, Ahmed and Akbar [31] numerically investigated the MHD Williamson nanofluid flow over an exponentially stretching surface.
Motivated by the above literature review, it has been observed that no study is so far reported which investigates the heat transfer characteristics of Williamson type fluid flow via an exponentially stretching curved surface in the presence of variable thermal conductivity. The prime aim of the present investigation is to study the heat flow of Williamson fluids over an exponentially stretching porous curved surface with a heat source. With the help of suitable similarity transformation, the governing PDEs are converted into ODEs. The highly non-linear ODEs are solved numerically by using Matlab code bvp4c. The impact of the curvature parameter, magnetic number, suction/injection parameter, permeability parameter, Prandtl factor, Eckert factor, non-linear radiation parameter, buoyancy parameter, temperature ratio parameter, Williamson fluid parameter, thermal conductivity parameter on velocity, pressure, and temperature profiles are observed by the plotted graphs. The variation in the skin friction and Nusselt number due to the involved physical parameters are observed through tables.

Problem Description
Here we have considered a steady, incompressible, laminar two-dimensional mixed convective Williamson type fluid flow under the impact of MHD to observe the attributes of the heat transfer mechanism along an exponentially stretching porous curved surface. Considering the coordinate system (r, s), where s is assumed along the fluid flow direction, and the r-axis is taken orthogonally to the flow, let R be the total radius of the curved surface. The physical scenario and the setup of coordinates are portrayed in Figure 1. The velocity field (stretching) of the surface is u = A 0 e s l , which is taken along the s-direction. An exponentially varying magnetic field B = B 0 e − s 2l is taken in the radial direction whereas, the temperature adjacent to the surface is T = T w . Under the above assumptions, the flow equations following [2,27] are given below: ∂p ∂r Here, the symbols u and v are in the s and r directions, respectively, representing the velocity components. g is the gravitational acceleration due to earth's attraction, β is the volumetric coefficient used for the thermal expansion, p is the dimensional pressure, and ρ, σ, µ are the fluid density, electrical conductivity of the fluid, and the dynamic viscosity, respectively. Γ is a positively taken time constant, and c p is the heat capacity of the particles, T is the temperature (fluid), q r is called the radiated heat flux, and k(T) is the variable thermal conductivity. Applying Rosseland approximation for radiation, we have: The T 4 Taylor series in terms of T ∞ is represented as: Using Equation (6) in Equation (5), we obtain: Here, σ represents the constant named Stefan Boltzmann and k represents the mean absorption parameter. The non-uniform source of heat or sink is represented by q and defined as: where D 1 , E 1 > 0 gives a rising push to the internal heat generation and D 1 , E 1 < 0 gives a rising push to internal heat absorption. k(T) is defined as: Here, is the variable thermal conductivity. Using Equations (7)- (9) in Equation (4), we obtain The accompanying boundary conditions are: We define the following similarity transformation to solve the governing equations: Here η, K, l, and θ w represent the dimensionless similarity variable, the curvature parameter, the reference length, and the temperature ratio parameter, respectively. f (η), θ(η), and P(η) represent the dimensionless velocity field, temperature distribution, and pressure profiles, respectively. Using the above transformations, in Equations (1)-(4), Equation (1) is satisfied identically, and the governing equation takes the following form: 4KP Fluids 2021, 6, 260 Equation (11) transfer into: by the elimination of P(η) from Equations (13) and (14): where is Williamson fluid parameter, A 0 e − s 2l is the Biot number. The quantities of physical importance, skin-friction (local wall-drag) coefficient C f , as well as the local Nusselt (heat flux) number Nu are defined by: The shear stress τ w and heat flux q w near the surface is defined as: Using Equation (12) in Equation (20), after simplification with Equation (19), we obtain the following: √ where Re s = u w l ν . The velocity field (stretching) of the surface is = , which is taken along the -direction. An exponentially varying magnetic field = is taken in the radial direction whereas, the temperature adjacent to the surface is = . Under the above assumptions, the flow equations following [2,27] are given below:

Result Analysis
We have considered a steady, incompressible, laminar, two-dimensional mixed convective Williamson type fluid flow under the direct influence of MHD to observe the properties of the heat transfer mechanism along an exponentially stretching porous curved surface. Equations (13), (15), (18), and (21), subject to the given sufficient boundary conditions (16) and (17), are solved numerically by using the MATLAB bvp4c code. The influence of dimensionless governing physical parameter of flow, i.e., Williamson fluid parameter λ, buoyancy parameter δ, Biot number Bi, permeability parameter Y, Eckert factor Ec, radiation parameter Rd, suction/injection parameter S, the magnetic number M, curvature K, Prandtl number Pr on velocity, pressure, and thermal profiles are shown through graphs and tables. The results are obtained by assigning K = 1.2; M = 0.5; = 0.5; S = 0.01; Y = 0.5; δ = 0.2; Ec = 1; Pr = 7; λ = 0.2; Bi = 0.1; θ w = 5; Rd = 0.1; D 1 = 0.02; D 2 = 0.01. For the validation of the code, a comparison has been made for the numerical values of the skin friction (local) coefficient for different values of the curvature parameter by fixing all parameters to zero with Kumar et al. [2] and Okechi et al. [27], and a good agreement has been found, as shown in Table 1. Numerical results/data of the wall drag and heat flux numbers for multiple numerical values of the dimensionless physical parameter are given in Tables 2 and 3. Table 2 shows that for elevated values of λ, a significant reduction in the skin friction (wall drag) is seen and the local Nusselt (heat flux) number declines slowly. This demonstrates that λ has strong influence on skin friction. For higher values of M, skin friction shows an augmented trend, whereas the local Nusselt factor decreases. Physically, the presence of the magnetic field lowers the motion of the fluid, which is applied perpendicular to the flow direction. For greater values of K, local skin friction (wall drag) and Nusselt (heat flux) number are on the downward trend. As we increase the permeability parameter, the local skin fraction increases whereas the local Nusselt number decreases. This effect is seen due to the rise in viscosity. In fact, the viscosity of a fluid weakens its flow velocity while boosting the local skin friction. The heat flux rate augments due to an increase in viscosity, and the surface temperature transfers their heat to the upper layers of the fluid. Subsequently, the wall temperature reduces. As we increase the suction/injection parameter, both the wall drag and the local Nusselt (heat flux) number show augmented behavior. By increasing the buoyancy parameter, the skin friction shows a certain decline, whereas the local Nusselt number shows enhancement.  Table 3 portrays the effects of the Biot number (Bi), Eckert factor (Ec), radiation parameter (Rd), Prandtl factor (Pr), non-uniformly induced heat sink, source parameters D 1 and E 1 , thermal ratio parameter (θ w ), thermal conductivity (variable) parameter ( ), on skin friction (wall drag) coefficient, and the local Nusselt (heat flux) number. Table 3 shows that an escalation in Ec results in the decline of the skin friction (wall drag) coefficient and local Nusselt (heat flux) number because Ec depends on the surface velocity. Ec augments mean the velocity field of the surface augments, and consequently, an anti-augmented trend in local skin friction (wall drag) and local Nusselt (heat flux) number is noted. By increasing the value of the Pr number, the skin friction and Nusselt number both show increasing behavior because the viscosity of the fluid increases and the flow rate decreases. Consequently, the decrease in the flow rate reports an augmented in-wall friction and wall thermal state. As we increase the thermal conductivity (variable) parameter , both the local Nusselt number and the local skin friction (wall drag) coefficient decline. The skin friction (wall drag) coefficient decreases for greater radiation parameter values and Biot numbers, whereas the Nusselt (heat flux) number increases. For augmented values of the heat source and sink parameters D 1 and E 1 both local skin friction (wall drag) coefficient as well as local Nusselt (heat flux) number decreases. Figure 2a,b shows the impact of Williamson type fluid parameter λ on the given velocity field f (η) and temperature θ(η) profiles. The significant influence of the Williamson fluid parameter imparted on velocity is noted in Figure 2a. It is noticed that for higher values of λ, the value of f (η) decreases, but Figure 3 shows that the rising numbers the of Williamson type fluid parameter λ yield an increase in the curve of the temperature θ(η). Physically, we can say that the greater values of λ mean more and more sufficient relaxation times, offering more resistance to the fluid flow. higher values of λ, the value of (η) decreases, but Figure 3 shows that the rising numbers the of Williamson type fluid parameter λ yield an increase in the curve of the temperature (η). Physically, we can say that the greater values of λ mean more and more sufficient relaxation times, offering more resistance to the fluid flow.  Figure 3a-c illustrate the effects of curvature K on f (η), P(η), and θ(η). It is evident from Figure 3a,b that the value of f (η) and P(η) increases with an enhancement in the value of the curvature parameter. Figure 3c portrays the effect of curvature K on the distribution of heat θ(η). By raising the value of the curvature factor, the graph shows a reduction in the value of heat distribution. As we increase K, the flatness of the stretching surface increases. Due to an increase in the flatness of the stretching surface, the flow velocity increases, and the temperature profile decreases because the resistance between the layers of the fluid reduces. Figure 4 shows the impact of Bi on heat distribution. As we increase the value of the Biot number, there is a certain rise in the temperature profile because the Biot number is a relation within the internal conductive resistance to the surface offered convective resistance. It means that as we increase the Biot number, the internal conductive resistance increases, making the cross-boundary layer thickness stronger. Figure 5 shows the impact of the Prandtl number given on the distribution of heat. We noticed from Figure 5 that we see a decrease in heat distribution for the rising values of the Prandtl number. The Prandtl number Pr is the ratio between the momentum and thermal diffusivity. Therefore, the Prandtl factor owes an inverse relationship with thermal diffusivity. So, by increasing the Prandtl number, there is a decrement in heat distribution. Figure 6a-c depicts the significant impact of magnetic field number M on velocity f (η), pressure P(η) and temperature profile θ(η). Figure 6a shows that with the augment in cumulative magnetic number M, the velocity f (η) depicts a reduction in values. Figure 6b,c show that an opposite trend is noticed. The increase in P(η) and temperature θ(η) result in an augmentation in the numerical values of M. This effect is seen due to increases in the Lorentz force. The Lorentz force acts perpendicular to the flow direction, which opposes velocity and increases pressure and temperature profiles. Figure 7 demonstrates the impact of θ w on the thermal profile. It is observed that when we raise θ w , the thermal profile also increases. Higher values of θ w means the temperature difference between T w and T ∞ increases. This temperature difference causes an increase in the temperature profile. Figure 8 illustrates the effect of non-linear thermal radiation factor Rd on the thermal profile. It is noted that when we increase radiation parameter Rd, the temperature profile also augments. Since the stretching surface is heated, the energy emitted in the form of electromagnetic radiation causes an increase in the temperature profile. sistance increases, making the cross-boundary layer thickness stronger. Figure 5 shows the impact of the Prandtl number given on the distribution of heat. We noticed from Figure 5 that we see a decrease in heat distribution for the rising values of the Prandtl number. The Prandtl number is the ratio between the momentum and thermal diffusivity. Therefore, the Prandtl factor owes an inverse relationship with thermal diffusivity. So, by increasing the Prandtl number, there is a decrement in heat distribution.   sistance increases, making the cross-boundary layer thickness stronger. Figure 5 shows the impact of the Prandtl number given on the distribution of heat. We noticed from Figure 5 that we see a decrease in heat distribution for the rising values of the Prandtl number. The Prandtl number is the ratio between the momentum and thermal diffusivity. Therefore, the Prandtl factor owes an inverse relationship with thermal diffusivity. So, by increasing the Prandtl number, there is a decrement in heat distribution.     Figure 10 portrays the impact of the Eckert number on heat convection phenomena. It is observed from the graph that by increasing the Eckert number, there is an increase in heat distribution. The Eckert factor is defined as the ratio between the kinetic energy as well as the boundary layer enthalpy. So, by augmenting the value of the Eckert number, there is a rise in the heat distribution value. Figure 11 shows the impact and significance of on the temperature profile θ(η). As we increase , the temperature profile also increases. Figure 12a-c give us the effect of the suction/injection parameter S on the velocity field, temperature distribution, and pressure profile. Figure 12a,b depict that as we increase the value of S, the velocity field, as well as temperature profile, declines. This effect is shown because by increasing the value of S, the suction/injection velocity increases, which is produced perpendicular to the direction of the flow, decreasing the motion of the fluid. The temperature profile decreases by increasing the value of S because the temperature near the surface is higher than the thermal state away from the surface by varying S; a reduction in temperature is shown in Figure 12b. Figure 12c depicts the pressure profile increasing as a result of the suction/injection factor on pressure from the interval (0,1); when the interval is (1,∞), the pressure profile decreases. Mathematically this effect is shown due to the impact of velocity because the pressure is the function of the square of the velocity (see Equation (13)).
( ), pressure ( ) and temperature profile ( ). Figure 6a shows that with the augment in cumulative magnetic number , the velocity ( ) depicts a reduction in values. Figure 6b,c show that an opposite trend is noticed. The increase in ( ) and temperature ( ) result in an augmentation in the numerical values of . This effect is seen due to increases in the Lorentz force. The Lorentz force acts perpendicular to the flow direction, which opposes velocity and increases pressure and temperature profiles.  Figure 7 demonstrates the impact of on the thermal profile. It is observed that when we raise , the thermal profile also increases. Higher values of means the temperature difference between and increases. This temperature difference causes an increase in the temperature profile. Figure 8 illustrates the effect of non-linear thermal radiation factor on the thermal profile. It is noted that when we increase radiation parameter , the temperature profile also augments. Since the stretching surface is heated, the energy emitted in the form of electromagnetic radiation causes an increase in the temperature profile.      Figure 10 portrays the impact of the Eckert number on heat convection phenomena. It is observed from the graph that by increasing the Eckert number, there is an increase in heat distribution. The Eckert factor is defined as the ratio between the kinetic energy as well as the boundary layer enthalpy. So, by augmenting the value of the Eckert number, there is a rise in the heat distribution value. Figure 11 shows the impact and significance of on the temperature profile ( ). As we increase , the temperature profile also increases.       Figure 11. Effect of on ( ). Figure 11. Effect of on θ(η).
Fluids 2021, 6, x FOR PEER REVIEW 14 of 18 Figure 12a-c give us the effect of the suction/injection parameter on the velocity field, temperature distribution, and pressure profile. Figure12a,b depict that as we increase the value of , the velocity field, as well as temperature profile, declines. This effect is shown because by increasing the value of , the suction/injection velocity increases, which is produced perpendicular to the direction of the flow, decreasing the motion of the fluid. The temperature profile decreases by increasing the value of because the temperature near the surface is higher than the thermal state away from the surface by varying ; a reduction in temperature is shown in Figure 12b. Figure 12c depicts the pressure profile increasing as a result of the suction/injection factor on pressure from the interval (0,1); when the interval is (1,∞), the pressure profile decreases. Mathematically this effect is shown due to the impact of velocity because the pressure is the function of the square of the velocity (see Equation (13)). Figure 13a-c reveal the effects of the permeability factor on velocity, temperature, and pressure profiles. As we increase the permeability factor , the velocity field decreases, and therefore, the temperature profile augments. The permeability factor is dependent on a viscosity, which means that increases the viscosity of the fluid under consideration , the flow of the velocity rate decreases, and the internal kinetic energy and heat transfer rates increase because of this effect. Figure 13c depicts the effect of the Figure 13a-c reveal the effects of the permeability factor on velocity, temperature, and pressure profiles. As we increase the permeability factor Y, the velocity field decreases, and therefore, the temperature profile augments. The permeability factor Y is dependent on a viscosity, which means that Y increases the viscosity of the fluid under consideration, the flow of the velocity rate decreases, and the internal kinetic energy and heat transfer rates increase because of this effect. Figure 13c depicts the effect of the permeability parameter on the pressure profile. The pressure profile increases from (0,1) and decreases from (1,5) for greater values of permeability parameter. Figure 14a-c depict the significant impact of the buoyancy parameter on the given velocity field, pressure profile, and thermal distribution. It is noticed from Figure  14a that by raising the value of , the velocity profile also categorically augments, whereas the temperature profile decreases (see Figure 14c). From a physical point of view, the buoyancy force acts on the curved surface and stimulates the movement of the particle near the surface; due to this effect, the velocity field augments, and the thermal state declines. An interesting phenomenon is shown in Figure 14b wherein a variation in from (0,1) pressure profile decreases and increases from (1,5). Figure 14a-c depict the significant impact of the buoyancy parameter δ on the given velocity field, pressure profile, and thermal distribution. It is noticed from Figure 14a that by raising the value of δ, the velocity profile also categorically augments, whereas the temperature profile decreases (see Figure 14c). From a physical point of view, the buoyancy force acts on the curved surface and stimulates the movement of the particle near the surface; due to this effect, the velocity field augments, and the thermal state declines. An interesting phenomenon is shown in Figure 14b wherein a variation in η from (0,1) pressure profile decreases and increases from (1,5).

Conclusions
In this article, heat transfer attributes of Williamson type fluid flows via an exponentially stretching porous curved surface are studied under the impact of thermal conductivity (variable) and convective conditions at the boundary. The governing partial differential equations are based on the Navier-Stokes theory. By applying suitable similarity transformation, governing problems are transformed into ordinary differential equations. After that, the solutions are obtained using the built-in MATLAB-based bvp4c package. The significant findings of this article are given below: • , , Υ, and reduce the velocity field. • , , Υ, and are the increasing function of pressure profile, however, it declines for . • The thermal distribution increases with increases in , , , , , , , Υ, and , whereas the thermal profile declines for greater values of , , , , and .

•
The skin friction (wall drag) coefficient shows a reduction for higher values of , , , , , , , , , and , whereas it increases for greater values of , Υ, , and .

Conclusions
In this article, heat transfer attributes of Williamson type fluid flows via an exponentially stretching porous curved surface are studied under the impact of thermal conductivity (variable) and convective conditions at the boundary. The governing partial differential equations are based on the Navier-Stokes theory. By applying suitable similarity transformation, governing problems are transformed into ordinary differential equations. After that, the solutions are obtained using the built-in MATLAB-based bvp4c package. The significant findings of this article are given below: