1. Introduction
Investigation of non-Newtonian fluidic models is a subject of curiosity for many researchers because of their nature, engineering, and industrial solicitations, including polymeric melt, crystal growing, dilute polymer solutions, drilling muds, cosmetic products, foods, coated sheets, and glass blowing. From a previous analysis, it is clear that many constitutional correlations have been put forward to flaunt the knotty traits of non-Newtonian fluids. Such entrenchment emerges from the view that the indispensable features of non-Newtonian liquids cannot be scrutinized using the classical Navier–Stokes equations, which are solely applicable to determine Newtonian fluid traits. Numerous researchers have investigated the behaviour of non-Newtonian fluids through various geometrical configurations. Among these geometries we have cylindrical surfaces, which play a significant role in polymer processing systems. Nazar et al. [
1] investigated mixed convective flow of micropolar fluids through a horizontal circular cylinder for the case of constant surface temperature. Madhavi et al. [
2,
3,
4] and Gaffar et al. [
5] scrutinized third-grade fluid flow through a horizontal cylinder surface. Gaffar et al. [
6] also explored boundary layer flow and the heat relocation of tangent hyperbolic fluid across a horizontal circular cylindrical body with first-order thermal and velocity slip. Nagaraj et al. [
7] investigated flow and heat relocation in non-Newtonian Eyring–Powell fluids over a horizontal circular cylindrical surface by taking into account suction/blowing and heat convective boundary constraints. Zokri et al. [
8] studied the flow of non-Newtonian Jeffrey fluids over a horizontal circular cylinder in the company of mixed convection and viscous heating effects.
The non-Newtonian fluid model for Williamson fluids was introduced by Williamson [
9], and it has been confirmed to be among the predominantly crucial non-Newtonian fluids due to its low viscosity and higher shear rate. In particular, the Williamson fluid model predicts that the effective viscosity will decline monotonically with an elevation in the shear rate and vice versa. This pseudo-plastic fluid deportment is one instance of the various ways in which this non-Newtonian liquid might be usable in modern technology and industry. The applicability of Williamson fluid in technologies and industries has generated a major effect on manufacturing processes. Thus, many scholars have been involved in the exploration of features of Williamson fluids in the improvement of the thermal features of this fluid. However, just like any other non-Newtonian fluid, the Williamson fluid is still incapable of fulfilling the demand for superior extremity heat transport due to its dissatisfying thermal conductivity. To upsurge thermal conductivity as well as enhance heat transport performance, the theory of adding nanoparticles (NPs) into the non-Newtonian fluid has proven to be one of the more contemporaneously innovative concepts in the recent past years. The mixture of NPs and the host fluid is called a nanofluid (NF) [
10]. Various models have been put forward to scrutinize NFs, of which the Buongiorno NF model [
11] is one of them. The Buongiorno [
11] NF model is a two-phase model that focuses on Brownian motion and thermophoresis impacts, and this model plays a significant role when the movement of NPs remarkably impacts the fluid flow. Another notable model is the Tiwari–Das NF model [
12], which is a single-phase model that is concerned with studying the volume fraction of NPs. Rashad et al. [
13] used the Buongiorno NF model to analyze the combined convection flow of a Newtonian NF across a horizontal circular cylinder immersed in a penetrable medium with heat-convective boundary conditions. Using the Tiwari–Das NF model, Zokri [
14] also studied the mixed convection flow of Jeffrey NFs past a horizontal cylinder by accounting for suction/injection and heat convective boundary constraints. Merkin [
15] initiated the study of fluid flow moving through a horizontal circular cylinder with mixed convection, and many studies were reported thereafter. Mixed convection plays a major role in enhancing the thermal features of heat transport. It also accounts for the general instance of convection that takes place in various industrial and technological solicitations including electronic device cooling, solar energy storage, drying technology, food processing, and float glass manufacturing.
In many engineering, technological, and industrial processes, the disparity between the surface and ambient (free stream) temperatures need not be small. In the instance of huge temperature variations, the fluid and heat transfer features are remarkably influenced. In such instances, accounting for linear Boussinesq approximation and classical linear Rosseland approximation in the fluid and heat transport equations, respectively, is inadequate, and that can lead to a reduction in the quality of the results. Consequently, it is mandatory to account for nonlinear Boussinesq approximations and nonlinear Rosseland thermal radiation in such instances so that accurate results are achieved. Practical examples of nonlinear thermal convection include engineering, geophysical, and astrophysical flows. On the other hand, nonlinear Rosseland heat radiation plays a fundamental role in manufacturing processes with larger temperature differences. EL-Zahar et al. [
16] employed the modified Tiwari–Das NF to investigate the magnetized mixed convective flow of hybrid NFs across a horizontal circular cylindrical surface by considering nonlinear thermal radiation. Basha et al. [
17] used the Buongiorno NF model to analyze the mixed convection flow of tangent hyperbolic NFs through a horizontal circular cylinder by accounting for nonlinear Boussinesq estimation. Activation energy (AE) is the minimal quantity of energy needed by chemical reactants to go through a chemical reaction. Mass transport procedures together with chemical reactions with Arrhenius AE emerge in areas of chemical engineering and food processing. Huang [
18] studied free convective flow via a horizontal circular cylindrical surface immersed in a penetrable medium with AE and suction/injection.
Bioconvection is a phenomenon whereby the flow and thermal features of the fluid are impacted by the existence of microbes. Bioconvection flows emerge in various fascinating solicitations in environmental systems, fuel cells, and biological polymer synthesis [
19]. The upward swimming microorganisms give rise to unstable density stratification, and because of that, bioconvection takes place in the fluid. In other words, bioconvection is compelled by instability that emerges due to self-propelled microorganisms swimming to and interacting with NPs and buoyancy [
20,
21,
22,
23]. The amalgamation of bioconvection and NF has yielded some exceptional results in micro-fluidic devices such as micro-channels and micro-rectors. Kuznestov and Avramenko [
24] initiated the use of NPs to scrutinize the bioconvection of gyrotactic microorganisms. Hill and Pedly [
25] provided a review study of the hydrodynamics of bioconvection. Through the use of Navier–Stokes equations, Allouit et al. [
26] investigated the bioconvection of gyrotactic microorganisms via a vertical cylinder. Mallikarjuna [
27] scrutinized the mixed bioconvective flow of a NF containing gyrotatic microorganisms from a vertical slender cylindrical surface. Rashad et al. [
28] investigated mixed the bioconvective flow of a Newtonian NF comprising gyrotactic microbes through a horizontal circular cylinder in the company of heat-convective and zero mass flux boundary conditions.
The literature review, particularly investigations by Rashad et al. [
13,
28], reveals that there is no study on the bioconvective flow of non-Newtonian Williamson NFs in thermal and mass relocation mechanisms over a cylindrical surface with the following research gaps that result in the novelties of the present analysis.
The complex rheology features of Williamson materials are examined using the modified Buongiorno NF model in the presence of an applied magnetic field, porous media, and variable fluid viscosity.
Nonlinear Boussinesq approximation and nonlinear Rosseland thermal radiation are incorporated in the momentum and energy equations, respectively, to cater for large temperature differences.
Heat generation/assimilation and variable thermal conductivity are incorporated in the heat equation for the characterization of energy transfer mechanisms.
The modified Arrhenius AE, chemical and microbial reactions, NP and microbial Brownian diffusions are considered in the NP concentration and motile microorganism equations for characterization of mass and motile microbes transfer mechanisms.
The main objective of the current analysis is to scrutinize the bioconvective flow, heat, mass, and motile microbes transfer features of the magnetized Williamson NF model containing motile gyrotactic microbes through a horizontal circular cylindrical surface encased in a porous medium. The consequences of nonlinear thermal convection, nonlinear heat radiation, heat generation/assimilation, varying fluid attributes, AE with chemical and microbial reactions, and Brownian motion for both NPs and motile microbes are further elaborated as novelties. There has been limited research conducted on the influence of microbial Brownian motion in conjunction with microbial reactions on fluid flow properties [
29,
30]. Moreover, the entire flow analysis has been subjected to first-order velocity slips, suction/injection, and heat-convective and zero mass flux constraints at the boundary. The consideration of mixed convective flow of non-Newtonian NF through a heated horizontal cylinder plays a fundamental role in geological, technological, and industrial processes including the thermal recuperation of oil, geothermal reservoirs, insulating wires, drag reduction, thrust bearing, assessment of aquifers, and underground nuclear waste storage. An efficient and accurate overlapping grid-based multi-domain bivariate spectral local linearization method (OMD-BSLLM) [
31,
32,
33] is executed to procure approximate solutions to the complex fluid flow model. The impact of various pertinent variables on the flow profiles and physical quantities of engineering attentiveness is scrutinized through graphical representations.
2. Mathematical Modeling
We contemplate the two-dimensional, laminar, steady, incompressible, and viscous nonlinear mixed convective boundary layer flow of Williamson NFs comprising gyrotactic microbes that swim from a horizontal circular cylindrical surface immersed in a porosity medium. The pictorial representation of the flow model and the physical coordinate system is depicted in
Figure 1, where
a is the radius of the cylinder. The
and
coordinates are, respectively, measured alongside the surface of the cylinder, resuming from the lower stagnation point, and normal to the cylinder. The uniform magnetic field intensity is enforced parallel to the fluid movement. Fluctuations in density for the buoyant expression are established by utilizing non-linear Boussinesq estimation, and the basic fluid and NPs are in a thermal equilibrium state. The surface of the cylinder is heated via convection from the heated fluid at the temperature
which produces a heat transport coefficient
The cylinder’s surface is maintained at a constant temperature
which is assumed to be greater than the free stream (ambient) temperature
for the warmed surface (aiding flow) and lower than the ambient temperature
for the cooled surface (reversing flow). The velocity of the external flow is represented by
where
is the ambient velocity. The surface of the cylinder is also maintained at a uniform density of motile microbes
, whereas the uniform NP volume fraction and density of motile microbes far-off from the cylindrical surface are signified by
and
respectively. Due to the absence of agglomeration and accumulation of NPs, the NF suspension is dilute (i.e homogeneous dispersion is attained). Furthermore, it is presumed that NPs have no impact on the direction and velocity of the gyrotactic microbes’ swimming. It is worth noting that the velocity for motile microbes, NPs, and the basic fluid is the same. For the characterization of heat transportation, features of nonlinear thermal radiation, heat sink/source, and variable thermal conductivity are involved in the heat (energy) equation. However, for the characterization of mass and motile microbes transportation, features of AE, chemical and microbial reactions, and NP and microbial Brownian motions are invoked in the equations of NP species concentration and motile microbe density conservation. In view of the aforementioned presuppositions, the entire fundamental equations of the flow problem take the form [
13,
28].
where
and
are velocity segments along the
and
axes,
and
represent the density of the fluid, NPs, base fluid, and microbes,
is the material constant,
signifies electrical conductivity,
stands for the intensity of the magnetic field,
stands for the dynamic viscosity that varies with temperature,
indicates thermal conductivity that varies with temperature,
g represents acceleration due to gravity,
is the penetrability of the porous medium,
and
signify linear and nonlinear volumetric thermal expansion coefficients,
is the average volume of the motile microbes,
stands for the heat sink/source coefficient,
represents the specific heat capacity,
signifies the NF heat capacity ratio,
is the thermophoretic diffusion coefficient,
is the Brownian motion factor,
is the coefficient of motile microbes diffusion,
and
are the chemical and microbial reaction rates,
is the Boltzmann constant,
p is the fitted rate constant,
is the coefficient of AE,
b is the chemotaxis constant, and
is the speed of the swimming cell. It is notable that the third term of the momentum boundary layer equation introduces the non-Newtonian fluid’s behaviour, while the fourth term represents contributions from the porous media and magnetic field. The remaining three terms on the right-hand side of the momentum equation account for buoyancy forces due to temperature variations, NP concentration, and motile microbes. In the heat equation, the second and third terms on the right-hand side capture the contributions from heat generation/absorption and thermal radiation flux based on the Rosseland approximation. However, the last term addresses the impact of nanoparticles on the thermal behavior of the fluid, encompassing phenomena such as thermophoresis and NP and microbial Brownian motions. The second term of the NP concentration and motile microbes equations incorporates the effects of AE, chemical reactions, and microbial reactions. With regard to the Rosseland diffusion approximation, the radiative flux is furnished as [
34,
35]
where
represents the Rosseland extinction coefficient,
signifies the Stefan–Boltzmann constant, and the term
is called radiative conductivity. Since features of first-order velocity slip, uniform suction/injection velocity, and thermal convective and zero mass flux conditions are considered at the boundary, the appropriate physical boundary conditions for the flow problem under consideration are given by
where
is the mass suction/injection velocity and
is the first-order velocity slip factor.
In view of Rashad et al. [
13,
28], the external flow velocity
for the boundary layer equations takes the form
where
is the free stream velocity. The next dimensionless variables are adopted in order to aid in the attainment of the numerical solutions [
28]
where
signifies the Reynolds number,
with
being the surface temperature excess ratio. The correctness of the NF flow and the rate of heat transmission can be achieved by considering the temperature-variant fluid viscosity and heat conductance. Consequently, the fluid velocity follows an exponential variation with temperature and the thermal conductivity fluctuates linearly with temperature, respectively, as follows [
36,
37,
38,
39]
where
is the viscosity variation parameter and
(
for fluids such as air and water,
for fluids such as lubricant oil) is the temperature-variant thermal conductivity parameter. Applying the non-dimensional variables into Equations (
1)–(5), we obtain the next non-similar PDEs
and the corresponding boundary conditions change to
The flow parameters appearing in Equations (10)–(14) are explicated as magnetic field parameter Williamson fluid parameter , Prandtl number thermal radiation parameter heat source/sink parameter porosity parameter velocity slip parameter thermal Biot number buoyancy ratio parameter bioconvection Rayleigh number nonlinear convection parameter NP Brownian motion parameter , microbial Brownian motion parameter thermophoresis parameter Schmidt number bioconvection Schmidt number AE parameter chemical reaction rate parameter bioconvection Peclet number bioconvection constant microbial reaction rate parameter and mixed convection parameter with being the Grashof number.
Equations (
9)–(14), along with boundary conditions (
15), can be expressed in dimensionless form by assuming that
and
[
13,
28], where
signifies the dimensionless stream function given by
and
Thus, the continuity equation is automatically satisfied and the momentum, energy, NP concentration, and microorganism equations convert into
with accompanying boundary conditions
where
is the suction parameter when
or the injection (blowing) parameter when
. It is noteworthy that near the lower stagnation point of the cylinder (i.e.,
),
and Equations (
16)–(19) reduce to the subsequent system of ordinary differential equations (ODEs):
and the boundary conditions diminish to
The significant physical quantities of engineering attentiveness are the local wall friction factor
, local Nusselt number
, and local density number of motile microorganisms
that can be derived using the transformations described above, and expressed in dimensionless form as follows:
where
(surface shear stress),
(surface heat flux) and
(surface mass flux).
3. Solution Method
Numerical solutions for the non-dimensional PDEs (
16)–(19) along with boundary conditions (
20) have been determined by applying OMD-BSLLM. This numerical method has been successfully implemented in solving dimensional PDEs for some boundary layer flow problems with different degrees of non-linearity and complexity [
31,
32,
33]. The choice of the numerical scheme for the present flow problem is justified by the following merits:
The numerical approach is remarkably accurate, convergent, stable, and resource-efficient when solving problems in both smaller and larger computational domains.
The numerical method uses fewer number of grid nodes and iterations to achieve results with good spectral accuracy.
This numerical method is well-structured, straightforward, and versatile to program in various computer software such as Matlab and Mathematica.
The scheme produces a less dense matrix system due to the use of the overlapping grid approach.
The numerical procedure for obtaining the numerical solutions to the present problem using the OMD-BSLLM is discussed in detail in the following section. The method is applied to linear equations, therefore it is necessary to begin by linearising the nonlinear equations using one of the linearisation techniques available in the literature. To this end, the algorithm called local linearisation method (LLM) is applied in the linearisation of the nonlinear terms in the dimensionless PDEs (
16)–(19). The LLM algorithm was identified by Motsa et al. [
40,
41,
42] as an efficient method for solving coupled systems of nonlinear ODEs and PDEs that model boundary layer flow problems. The LLM algorithm used in the simplification and decoupling of Equations (
16)–(19) is summarized in
Appendix A. The resulting linear iterative scheme is as follows:
where
represent current and previous iterations, and the variable coefficients have been obtained using the LLM approach (see
Appendix B). The next step involves discussing the decomposition of the time domain
which is the solution domain for the dimensionless PDEs. The time interval
is then segmented into
non-overlapping sub-domains indicated as follows:
where each sub-domain is further partitioned into
Chebyshev collocation points (see
Figure 2). The PDEs are solved in the
non-overlapping sub-domains
by starting with the easily calculated initial solutions in the first sub-interval
. These initial solutions, which are initial conditions for the first interval, are obtained by solving the system of ODEs (
21)–(
24). For the remaining sub-domains, we use the solutions produced at the last grid point of each sub-interval as initial conditions for the next sub-interval. This means that the solution at the last grid point of
is used as an initial condition in
In accordance with the validity of the spectral collocation method (SCM), which is in the interval
the decomposed domain
has to be transformed to
using the following linear transformation:
Because of this transformation, we now have a new collocation variable
which is defined as [
43]
The complete grid in the time variable (
) takes the following order:
Next, we discuss domain decomposition for the spatial variable
Unlike the variable
which is defined in the finite domain
the space variable
is specified in the semi-finite domain
In order to apply the method, this semi-finite domain has to be truncated to a finite domain
where
is a finite value that is handpicked with the intention that it conforms to the far field boundary constraints. The resulting truncated domain is then divided into
overlapping sub-domains indicated as
where each sub-domain is subsequently partitioned into
Chebyshev collocation points (see
Figure 3). Yang et al. [
44] reported that using a different number of collocation points in the various sub-intervals can lead to difficulties in the implementation of the method. In addition to that condition, the length for each sub-domain must be kept the same for all sub-domains, especially if a linear transformation is utilized. This length is given by
(see
Appendix C). Again, the decomposed domain
has to be transformed to
as well, using the following linear mapping:
and because of this transformation, we have a new collocation variable
which is defined as
The complete grid in the space variable (
) takes the following order:
The numerical solutions at different time intervals can be differentiated using
, and
We remark that for the time variable, equations are solved independently in each sub-domain, but for the space variable, equations are solved contemporaneously across the whole domain of integration. Thus, the linearised iterative scheme (
27)–(30) takes the form
In the first interval
the initial conditions
, and
are used to obtained numerical solutions
, and
Then, for the remaining time sub-intervals
, the continuity conditions are given as
are used to generate numerical solutions
, and
in the other sub-intervals. The procedure for generating the numerical solution is carried out by approximating the desired solutions, such as
, using the bivariate Lagrange interpolating polynomial that appears as
For each space sub-domain
, the first, second, and
order derivative matrices in the spatial direction are evaluated at the collocation points
for
as follows [
33]:
where the Chebyshev differential matrix
and
in the standard Chebyshev differential matrix of the first order with dimension
However, the matrix
has dimension
where
is the combination of all collocation points used in the whole spatial region. On the other hand,
is the matrix-vector function with dimension
We note that the
derivative matrix is obtained via matrix multiplication. For each time sub-interval
, the first-order derivative matrix in the time direction is evaluated at the collocation points
for
as [
33]
where
and
in the standard Chebyshev differential matrix of the first order with dimension
It is noteworthy that the matrix-vector function
has dimension
, since it accounts for the solution in the entire spatial domain. A similar procedure is utilized in the approximation of the other unknown functions
, and their corresponding derivatives at the collocation points
. Because of the execution of the overlapping grid strategy, the Chebyshev differential matrix
contains many zero elements and becomes less dense, and the structure of the matrix is illustrated in [
32,
33,
38,
39] (see
Appendix D). The crucial idea behind the use of SCM lies in the approximation of continuous derivatives with discrete derivatives. Consequently, by replacing the continuous derivatives in Equations (
39)–(42) with the discrete derivatives discussed in the previous paragraphs and making use of the initial conditions, we obtain
where the right-hand side of the equations are explicated as
and
The matching boundary constraints are also evaluated at the collocation points, and they take the form
In the form of an
matrix arrangement, Equation (
48) takes the form
where
and
is an
unit matrix. The additional equations will also give a similar matrix arrangement of the same dimension. The resulting discrete boundary conditions are incorporated into the main diagonal sub-blocks of each matrix system to yield a new system of linear algebraic equations. The later matrix system is then solved iteratively, resuming with initial approximations that are picked in such a way that they conform to the boundary conditions. The compact form of matrix Equations (
48)–(51) is given by
The functions used as initial guesses are given by
4. Results and Discussion
A comprehensive parametric study is carried out and numerical results are reported in this section through graphs. The numerical computation is performed for each parameter by assigning fixed values to other related parameters as
[
13,
28],
[
28],
[
28],
[
13,
28],
[
13,
28],
[
13,
28],
[
28],
[
28],
[
13,
28],
[
28],
[
4],
[
18],
[
13,
28],
[
4,
17],
[
18]
[
17],
[
16],
[
17],
[
28],
[
37],
[
37],
[
39],
[
18],
[
18],
[
17],
[
13,
28], and
[
29]. To corroborate the correctness of the present numerical method, the results for the heat transfer rate
were contrasted against those disclosed by Merkin [
15], Nazar et al. [
1], Rashad et al. [
13], Rashad et al. [
28], and Zahar et al. [
16] in
Table 1.
Table 1 demonstrates that the numerical values in all considered cases were in substantial coherence. Hence, the selection of the present overlapping grid-based numerical scheme can be rationalized.
The convergence of the OMD-BSLLM is tracked by monitoring the solution errors, as defined in [
32]. These errors represent solution-based discrepancies and quantify the number of accurate digits in the approximate solutions at the
r-th iteration level. If the numerical scheme is converging, it is anticipated that the error norms will decrease as the number of iterations increases.
Figure 4 illustrates the evolution of the solution errors for the approximate numerical solutions of
and
over iterations in both MD-BSLLM and OMD-BSLLM. The consistent reduction in all solution errors suggests that the numerical methods converge. Full convergence is considered attained when the convergence plots start to level off. It is evident that complete convergence is reached after approximately five iterations for all solutions in the MD-BSLLM approach and six iterations for all solutions in the OMD-BSLLM technique, with solution errors approaching
and
in the MD-BSLLM scheme and OMD-BSLLM technique, respectively. The minimal errors observed with the OMD-BSLLM approach validate its superior accuracy compared to the MD-BSLLM scheme. The precision of the overlapping grid-based scheme (OMD-BSLLM) can be assessed by evaluating the residual errors, as defined in [
32]. These errors quantify how closely the numerical solutions approximate the true solution of the PDEs (
16)–(19).
Figure 5 is plotted to demonstrate residual error graphs at different time levels for both OMD-BSLLM and MD-BSLLM (no overlapping) schemes. In all the graphs, it is clear that each graph displays sufficiently small residual errors, and that is an indication of good accuracy in both numerical schemes. The accuracy is also noted to be uniform throughout the time domain, and that is one of the benefits of using the multi-domain approach in the implementation of SCM. Nevertheless, the OMD-BSLLM scheme proves to be more accurate than the MD-BSLLM approach in all the graphs. Furthermore, using more overlapping sub-domains with a smaller number of spatial collocation points yields more precision in the OMD-BSLLM approach.
The stability of the OMD-BSLLM scheme is assessed by analyzing the condition numbers of the coefficient matrices associated with the system of linear algebraic equations that need to be solved. The condition number of the matrix gauges the degree of sensitivity of the solution to variations in the input data and to the round-off errors incurred during the solution process. The condition numbers for the coefficient matrices
are presented in
Table 2 when the number of collocation points and overlapping sub-intervals are varied. The overlapping grid-based numerical scheme is noted to yield small condition numbers compared to the non-overlapping grid-based scheme. The use of the overlapping grid approach results in coefficient matrices that are better conditioned compared to those generated by the non-overlapping grid-based calculations. The small condition numbers observed with the overlapping grid-based scheme indicate that the problem is well-conditioned, suggesting that the matrix representing the problem is not close to being singular. This characteristic is desirable as it ensures the stability, accuracy, and reliability of the numerical solution to the problem. As the number of overlapping sub-intervals increases and the number of space collocation points decreases, the condition numbers also decrease. This trend suggests that maximizing the number of overlapping sub-intervals while reducing the number of space collocation points further enhances the stability and accuracy of the numerical solutions. These findings affirm that the OMD-BSLLM scheme is the preferred numerical method for solving fluid flow problems akin to the one addressed in the current study.
For brevity, the consequences of other flow parameters such as
and
on the fluid flow attributes and design quantities of engineering attentiveness have been omitted in the current study. This is because the contributions of these flow parameters are similar to the ones that were adequately discussed in Rashid et al. [
13,
28] in the scenario of Newtonian fluids and the absence of variable fluid properties, nonlinear heat convection, magnetic field, heat sink/source, nonlinear thermal radiation, bioconvection Brownian diffusion, activation energy, chemical and microbial reactions, suction/injection, and first-order velocity slips. The current study seeks to extend the works of Rashid et al. [
13,
28] by analyzing, from a mathematical point of view, the fluid flow properties, heat, mass, and motile microbe transfer phenomena when the above features are incorporated in the present bioconvection flow model.
Figure 6,
Figure 7,
Figure 8,
Figure 9,
Figure 10,
Figure 11,
Figure 12,
Figure 13,
Figure 14,
Figure 15,
Figure 16,
Figure 17,
Figure 18,
Figure 19 and
Figure 20 illustrate the impact of various parameters stemming from the aforementioned novel features, including the Williamson fluid
, nonlinear convection
, suction
, velocity slip
, variable viscosity
, magnetic field
, thermal radiation
, temperature ratio
, variable thermal conductivity
, heat source
, chemical reaction
, AE
, microbial reaction
, microbial Brownian motion
, and motile microbe
parameters, on dimensionless velocity
, temperature
, NP concentration
, density of motile microbes
, wall friction factor
, Nusselt number
, and density number of motile microbes
. These pertinent parameters are chosen in the ranges:
= 0, 1, 2, 3;
= 0.1, 0.5, 1, 1.5;
= 0.1, 0.3, 0.5, 0.8;
= 0.1, 0.3, 0.5, 0.8;
= 0.1, 0.5, 1, 1.5;
M = 0.1, 1, 5, 10;
= 0.3, 0.6, 1, 1.5;
= 1.1, 1.3, 1.5, 1.8; 0.1, 0.5, 1, 1.5;
Q = 0.1, 0.3, 0.5, 0.8;
= 0.1, 1, 2, 3, 5;
E = 0, 1, 2, 3;
= 0.1, 1, 3, 5;
= 0.2, 0.3, 0.5, 0.8; and
= 0.1, 0.3, 0.5, 0.8. The selected range of values for the present model align with those commonly found in a typical nanofluid. For parameters exhibiting a broad spectrum, we opted for values that ensured the stability of the numerical scheme.
Figure 6 is plotted to demonstrate the impact of the Williamson fluid parameter
and nonlinear thermal convection parameter
on the dimensionless velocity. From
Figure 6a, it is clear that velocity dissemination is diminished with increasing values of
The variable
signifies the ratio of Williamson fluid relaxation time to retardation time. When the Williamson fluid parameter gets larger, the relaxation time of the fluid is noted to be lifted. As a result, resistance in the fluid flow direction is produced, which in turn reduces the NF velocity. Basically, because of the high relaxation time, fluid particles take more time to migrate to their original position; thus, fluid viscosity is improved and the velocity of the fluid particles is diminished. Conversely, velocity dissemination is accelerated with mounting values of
. For larger
, the disparity between the surface and free stream (ambient) temperatures
is increased, instigating a rise in the fluid velocity.
The consequence of the suction parameter
and first-order velocity slip parameter
on the velocity dissemination is portrayed in
Figure 7.
Figure 7a demonstrates that the velocity boundary layer’s thickness is diminished with the mounting values of the suction parameter. This deportment is owed to the fact that suction upshot plays a fundamental role in the removal of the fluid from the flow system; thus, momentum boundary layer thickness is remarkably decelerated.
Figure 7b shows that dimensionless velocity is maximized with the growing velocity slip variable. This is on account of the acceleration in the fluid flow and the corresponding momentum boundary layer thickness as the velocity slip effect is increased. The acceleration of the fluid velocity component with improved velocity slip impact was also disclosed by Gaffar et al. [
6] and Prasad et al. [
45] in the absence of microbes.
The variation of dimensionless velocity against the transverse coordinate
when the variable viscosity parameter
and the applied magnetic field parameter
are varied is illustrated in
Figure 8. It is clear that the velocity field is improved with rising values of the variable viscosity parameter (see
Figure 8a). Physically, increasing values of
tend to augment the temperature variation
consequently, the Williamson fluid bond is weakened and the intensity of the dynamic viscosity of the non-Newtonian material is minimized. As a result, the momentum boundary layer thickness is improved, and the fluid velocity is also accelerated. On the other hand, the dimensionless velocity of our nonlinear mixed convective flow is seen to be increased with mounting values of applied magnetic field parameter. This behaviour is opposite to that of the free convective flows (see Gaffar et al. [
6]), where the velocity is declined by rising values of
From a physics perspective, in an electrically conducting fluid, the magnetic field operating in a transverse direction to the geometry generates a resistive force known as the Lorentz force. The resulting Lorentz force inhibits the velocity of the fluid, as is evident in the free convective flows. In the present problem, it is noteworthy that the velocity component portrays the mounting behaviour because of the existence of nonlinear mixed convection. The similar behaviour is adequately explained by Basha et al. [
17] in the absence of gyrotactic microorganisms.
Figure 9 presents the effect of the thermal radiation parameter
and surface temperature excess ratio
on the dimensionless temperature curves. From the first figure, it is noted that the temperature curve is enhanced by rising values of the thermal radiation parameter. This is because increments in thermal radiation energize the boundary layer and upsurge the heat radiation energy; consequently, both the fluid temperature and thermal boundary layer width are also improved. The Rosseland radiation absorption is enhanced by the rising values of the heat radiation parameter
; thus, radiative heat transmission rate to the fluid is upsurged, leading to an improvement in the temperature dispersion. On the other hand, the non-dimensional temperature is also increased by the escalating values of the temperature ratio parameter since the parameter
improves the fluid’s thermal state. From a physics point of view, the greater values of the surface temperature excess ratio imply an upsurge in the temperature difference
consequently, the temperature of the fluid is elevated.
Figure 10 illustrates the ramifications of the variable thermal conductance parameter
and heat generation parameter
on the temperature curves.
Figure 10a reveals that the temperature curve improves with the escalation in
which is owed to the excess amount of heat transmitted from the surface of the cylinder to the fluid. The parameter
is directly proportional to the temperature difference. As a result, an elevation in the variable thermal conductance parameter causes a remarkable potential for thermal energy to convey more heat energy to the NPs, instigating a rise in the fluid temperature.
Figure 10b shows that temperature distribution mounts with the rising values of the heat source parameter. It is worth pointing out that
represents the amount of heat produced within the system during the fluid flow. Thus, an increment in
Q means that more heat is absorbed by the fluid, leading to an enhancement in the thermal boundary layer’s thickness and fluid temperature.
Figure 11 reveals the impact of the non-dimensional AE parameter
and chemical reaction rate parameter
on the non-dimensional NP species concentration distributions. From
Figure 11a, it is clear that increasing the non-dimensional chemical reaction rate lessens the thickness of the non-dimensional concentration dissemination. The thinning in the non-dimensional species concentration dissemination is continuously occurring in conjunction with an enlarged non-dimensional concentration gradient. On the other hand, the rise in the non-dimensional AE impact decreases the non-dimensional wall concentration gradient. The physical justification for such a behaviour is that high AE slows down the chemical reaction process. As a result, both the non-dimensional species concentration profiles (see
Figure 11b) and the thickness of the concentration boundary layer improve, leading to a reduction in the non-dimensional wall concentration gradient.
The impact of the AE parameter
and microbial reaction rate parameter
on the non-dimensional density of the motile microbes is portrayed in
Figure 12. It is perceived from
Figure 12a that the density of the motile microbes is reinforced by the rising values of the AE parameter. Since the microbial reaction process is improved by using AE, the motile microbe concentration profile is augmented to a greater extent. On the other hand,
Figure 12b illustrates that the density of motile microbes declines with the augmenting values of the microbial reaction rate parameter. With an elevation in
the microorganisms’ diffusivity shrinks, for which less motile microorganism transfer occurs. In this phenomenon, the concentration of motile microorganisms drops as the values of
increase.
Figure 13 reveals the non-dimensional microorganism concentration dissemination for the different values of the motile microbe parameter
and microbial Brownian motion parameter
. It is worth pointing out that the microbial Brownian motion is implied by the unrestricted mobility and unpredictable behaviour of motile microbes. The decline in the gyrotactic microorganism concentration profile for high
indicates fewer motile microbes for the reaction (see
Figure 13a). From a physics point of view, when
gets larger, the random collision between motile microbes increases, which leads to a reduction in motile microorganism transfers. As a result, an elevation in
resembles a decrement in the motile microbe concentration characteristics. In the other figure, there is an enhancement in the magnitude of motile microbes, in particular near the wall of the cylinder, when the bioconvection constant is raised. From formula
we note that the motile microbe constant correlates the ambient density of the motile microbes to the density difference across the boundary layer. When
upsurges, there is a larger density gradient across the boundary layer regime, which promotes the propulsion of motile microbes from the surface to the bulk flow. This results in an enhancement in the density of the motile microbes, in particular near the wall of the cylinder.
The ramification of the Williamson fluid parameter
and nonlinear thermal convection coefficient
on the wall friction factor is disclosed in
Figure 14.
Figure 14a demonstrates that the surface drag coefficient is a decreasing function of the changing Williamson fluid parameter. Larger values of
imply that a long relaxation time is needed to yield excess fluid movement resistance, thus increasing the skin friction factor.
Figure 14b shows that wall shear stress is also improved with a rise in the nonlinear temperature parameter since an upsurge in
raises the density of the fluid flow, which in turn augments resistance to the flow.
Figure 15 illustrates the upshots of the suction parameter
and velocity slip parameter
on the surface drag coefficient. It is perceived from
Figure 15a that a greater suction parameter leads to a decrement in the skin friction factor. This demonstrates that an improvement in the cylinder surface’s porosity causes an increment in the fluid flow resistance. In the other figure, surface shear stress is also diminished by increasing the values of
, since the fluid flow is decelerated along the cylinder when
is raised. As noted in Prasad et al. [
45], there is progressive movement in the peak shear stress locations far away from the lower stagnation when
is elevated. Thus, the repercussion of the wall slip is noted to be significant on the boundary layer features of the Williamson fluid flow from the horizontal cylinder.
Figure 16 highlights the consequences of the variable viscosity parameter
and applied magnetic field parameter
on the wall friction factor. From the first figure, it is clear that the skin friction factor is reduced by the rise in
because of a lower resistance to flow on the cylinder surface. When there is less resistance to flow, the fluid moves more freely along the surface, causing a reduction in the shear stress and skin friction coefficient. On the other hand, a rise in the magnetic field parameter augments the surface drag coefficient. This is because the Lorentz force monitors the rate at which the fluid particles tend to move; thus, the drag coefficient increases at the level of the cylinder’s surface.
The upshots of the heat radiation parameter
and temperature ratio parameter
on the dimensionless Nusselt number are depicted in
Figure 17. From
Figure 17a, it is clear that the elevation in
suggests a prominent increment in the local Nusselt number. This deportment is in harmony with the fact that heat transportation becomes superior with the radiation effect.
Figure 17b illustrates that the heat transfer rate escalates with increment in the surface temperature excess ratio since the presence of the temperature ratio parameter
contributes towards enhancement of the temperature of the fluid.
Figure 18 illustrates the distribution of the local Nusselt number for various variable thermal conductance
and heat generation parameters
. The local Nusselt number is enhanced by rising values of
as seen in
Figure 18a. This behaviour is due to the fact that the cylinder distributes additional heat by intensifying
which in turn accelerate the heat distribution in the system as well as the rate of heat transport.
Figure 18b demonstrates that escalation in
Q gives rise to higher Nusselt numbers. This is because the mechanisms of heat source produce a layer of warmed fluid near the cylindrical surface; thus, the temperature of the fluid surpasses the temperature of the cylinder’s surface.
The variation in the dimensionless local density number of the motile microbes with the AE parameter
and microbial reaction parameter
is depicted in
Figure 19. It is evident that the density number of motile microbes drops when the dimensionless AE parameter is increased. This signifies that as the Reynolds number escalates, the diffusion rate prevails over the motile microbe transfer rate. The lower the density number of motile microbes, the more AE there is. Higher AE implies that more energy is needed for the motile microbes to have a fruitful collision. On the other hand, the density number of motile microorganisms is noted to be remarkably enhanced by strong motile microbial reactions.
Figure 20 elucidates the impact of the bioconvection constant
and microbial Brownian motion parameter
on the density number of motile microbes.
Figure 20a reveals that as
augments, the density number of motile microbes is enhanced. This might be attributed to the fact that motile microbes in the Williamson fluid are moving haphazardly. The other figure demonstrates that the motile microbe transfer rate is reduced with the growing values of the motile microbes parameter.