Finite Difference Method to Evaluate the Characteristics of Optically Dense Gray Nanofluid Heat Transfer around the Surface of a Sphere and in the Plume Region

Abstract

The current research study is focusing on the investigation of the physical effects of thermal radiation on heat and mass transfer of a nanofluid located around a sphere. The configuration is investigated by solving the partial differential equations governing the phenomenon. By using suitable non-dimensional variables, the governing set of partial differential equations is transformed into a dimensionless form. For numerical simulation, the attained set of dimensionless partial differential equations is discretized by using the finite difference method. The effects of the governing parameters, such as the Brownian motion parameter, the thermophoresis parameter, the radiation parameter, the Prandtl number, and the Schmidt number on the velocity field, temperature distribution, and mass concentration, are presented graphically. Moreover, the impacts of these physical parameters on the skin friction coefficient, the Nusselt number, and the Sherwood number are displayed in the form of tables. Numerical outcomes reflect that the effects of the radiation parameter, thermophoresis parameter, and the Brownian motion parameter intensify the profiles of velocity, temperature, and concentration at different circumferential positions on the sphere.

1. Introduction

Nanotechnology is a key component of the industrial technological revolution of the twenty-first century, and it has become the interest of several researchers and engineers. Several industries, including petrochemical, textile, chemical, food processing, air conditioning, renewable energies, microelectronic cooling, etc., are facing the problem of high energy consumption in the thermal processes. Due to their enhanced thermophysical properties, the use of nanofluids is one of the non-destructive solutions that allows for the improvement in heat transfer performances. When nano-sized particles are suspended into a base fluid, the engineered suspension is termed ‘nanofluid’. The nanoparticles used to prepare the nanofluid are commonly made up of carbon nanotubes, oxides, metals, and carboids. The commonly used base fluids are ethylene glycol, water, and oil.

By dispersing metallic nanoparticles in traditional heat transfer fluids, Choi and Eastman [1] invented a new class of heat transfer fluids. The resulting nanofluids offered the opportunity of improving the heat transfer due to their enhanced thermal conductivities compared to conventional fluids. Brownian diffusion and thermophoresis are significant slip mechanisms in nanofluids. Buongiorno [2] created a two-component, four-equation, nonhomogeneous equilibrium model for heat, momentum, and mass transfer in nanofluids based on this discovery. The unexpected increases in the heat transfer coefficient cannot be explained by nanoparticle dispersion, according to a nondimensional examination of the equations. An electrically conducting incompressible nanofluid fluid flowing across a porous rotating disc in the presence of an externally imposed uniform vertical magnetic field was the subject of the study of Rashidi et al. [3]. This research has potential for nuclear propulsion space vehicles as well as spinning magnetohydrodynamic (MHD) energy producers for novel space systems and thermal conversion techniques. The concept of a three-dimensional flow of viscous nanofluid in the presence of partial slip and thermal radiation effects was introduced by Hayat et al. [4]. The flow is caused by a stretching porous surface. Alumina was chosen as the nanoparticle and water as the base fluid. In their research, fluid conducts electricity when a magnetic field is applied. The effects of non-linear heat radiation, slip conditions, and sheet stretching linearly on the phenomena of three-dimensional MHD nanofluid flow have been investigated. The radiative heat transfer of a nanofluid flow caused by the stretching and shrinking of a surface has been considered by Pal et al. [5]. The authors also studied the influences of dissipation and magnetics. Bioconvection heat transfer and radiative magneto-hydrodynamic nanofluid flow with a homogenous, heterogeneous chemical reaction past paraboloid geometry was the focus of the work of Makinde et al. [6]. Ahmad and Nadeem [7] studied the mechanism of hybrid nanofluid stagnation point flow produced by a rotating disk while considering the impact of an applied magnetic field. The discussion of the phenomenon of the natural convective flow of nanofluids in a plume-region was documented by Ashraf et al. [8]. Khan et al. [9] examined the effects of heat generation and applied magnetic field on nanofluid-free convection at several circumferential points of the sphere and in the plume-region in [9]. Effects of Hall currents, activation energy and chemical reaction around a rotating cone are taken into discussion by Ahmad and Nadeem [10].

The convective heat transfer combined with thermal radiation effects has been and is still the interest of the research community due to its various applications in industry and engineering, such as boilers, rocket engines, the cooling of the nuclear reactors, and thermal insulation. Hossain and Thakkar [11] explored the model of convective flow over the vertical plate kept at uniform temperature by considering radiation influence. Damseh [12] performed a numerical analysis on the magnetohydrodynamic-mixed convective flow over an isothermal surface immersed in a porous medium. Zehamatkesh [13] presented a study on the thermal radiation’s impact on free convective flow through a permeable enclosure. Mondal et al. [14] investigated the radiative convective flow across a vertical plate immersed in a porous media. Muhammad et al. [15] considered micropolar fluid flows past a porous moving surface situated in a spongy medium while being affected by heat generation, magnetic fields, radiation, and chemical reactions. Bhattacharyya et al. [16] examined the micropolar fluid flow with radiation effects over a stretching sheet. Pal et al. [17] were interested in the micropolar fluid’s stagnation point flow when radiation and a magnetic field were present. Moradi et al. [18] examined the mixed convection flow down the sloped surface with the influence of radiation. Sheikholeslami et al. [19] paid their attention to the energy transport mechanism, in nanofluid-filled semi-annulus in the presence of a magnetic field. Hussain et al. [20] reported the effects of radiation and transpiration on the boundary layer flow of a micropolar fluid past stretching sheet. The boundary layer flow over an exponentially stretched sheet with transpiration, radiation, and slip condition effects were explored by Mukhopadhyay et al. [21]. The thermally stratified radiative flow of third-grade fluid flow over a stretching surface was examined by Hayat et al. [22]. Parkesy et al. [23] studied the influence of radiation on transient magnetohydrodynamic micropolar fluid flow between porous vertical channels with imposed third-kind boundary conditions. According to Uddin et al. [24], the magnetohydrodynamic heat transport in a micropolar fluid past a wedge with a hall and ion slip current is affected by thermal radiation, heat generation, and absorption. The stagnation point flow of nanofluid with magnetohydrodynamic and radiation effects over a stretching sheet with velocity and thermal slip influences were proposed by Haq et al. [25]. Mohammad et al. [26] considered the mixed convection flow of micropolar fluid caused by an unsteady stretching sheet with viscous dissipation effects. A magnetohydrodynamic non-Darcian convective flow over a stretching sheet in micropolar fluid under the impact of solar rays was discussed by Mohammd et al. [27]. The discussion on radiation effects combined with variable density on magnetohydrodynamic Sakiadis flow was conducted by Abbas et al. [28]. Ashraf et al. [29] studied numerically the periodic-mixed convective flow along a sphere at different points of the surface. Ahmad et al. [30] studied the effects of the catalytic chemical reaction modification on the mixed convective flow along a curved surface. Spherical geometries exist in several industrial and engineering applications. By considering various fluid parameters, thermophoretic and Brownian motion effects, convective heat transport along the surface of a sphere and other geometries has been thoroughly anticipated in Refs. [31,32,33,34,35,36,37,38,39,40,41,42,43].

From the above-mentioned introduction and literature survey, it is shown that an enormous amount of research on nanofluids using the single-phase model and two-phase model have been documented due to the practical applications of nanofluids. Different geometries along with diverse physical characteristics of the fluids were considered to study the physical behavior of nanofluids, as in Refs. [1,2,3,4]. However, from the best of our knowledge, the study focusing on the two-phase model (Buongiorno model) through a sphere and in the created plume-region above the sphere due to the eruption of fluid has not been discussed before the present work. In the current work, linear thermal radiation effects of optically dense gray nanofluid are discussed. To obtain accurate solutions, an efficient numerical technique based on finite difference method is efficiently utilized. The outcomes of the pertinent parameters appeared in the model are portrayed in graphs and tables and are discussed with physical interpretations.

2. Flow Analysis for Region-I

The studied configuration consists of a two-dimensional, steady, viscous, and incompressible nanofluid flow around the sphere and in the plume-region. The plume-region is established just above the sphere due to eruption fluid. The thermal radiation effects along with thermophoresis and Brownian motion effects are considered in the following model. The base fluid is water and due to radiative heat flux, it is considered to be optically dense gray. The surface of the sphere is maintained at the temperature ${\widehat{T}}_w$ that is higher than the ambient temperature ${\widehat{T}}_{\infty }$. The mass concentration on the surface of the sphere is ${\widehat{C}}_w$ and this concentration is higher than the concentration in the free stream region ${\widehat{C}}_{\infty }$. The axes $\widehat{x}$ and $\widehat{y}$ are taken toward and normal to the direction of flow, respectively, as presented in Figure 1. The velocity components are designated by $\widehat{u}$, and $\widehat{v}$, respectively. Three regions make up the flow’s schematic diagram: region-I is the spherical boundary layer, region-II is the eruption region, and region-III is the plume region. Following [1,2,3,4,8,9], the governing equations are as follows:

Law of conservation of mass (Equation of continuity),

$$\frac{\partial \left(\widehat{r}\widehat{u}\right)}{\partial \widehat{x}}+\frac{\partial \left(\widehat{r}\widehat{v}\right)}{\partial \widehat{y}}=0.$$

(1)

Boundary layer form of momentum equation

$$\widehat{u}\frac{\partial \widehat{u}}{\partial \widehat{x}}+\widehat{v}\frac{\partial \widehat{u}}{\partial \widehat{y}}=\nu \frac{{\partial }^2\widehat{u}}{{\partial }^2\widehat{y}}+g{\beta }_T\left(\widehat{T}-{\widehat{T}}_{\infty }\right)Sin\frac{\widehat{x}}{a}+g{\beta }_C\left(\widehat{C}-{\widehat{C}}_{\infty }\right)Sin\frac{\widehat{x}}{a},$$

(2)

Boundary layer form of energy equation

$$\widehat{u}\frac{\partial \widehat{T}}{\partial \widehat{x}}+\widehat{v}\frac{\partial \widehat{T}}{\partial \widehat{y}}={\alpha }_m\frac{{\partial }^2\widehat{T}}{{\partial }^2\widehat{y}}+\tau \left\{D_B\frac{\partial \widehat{C}}{\partial \widehat{y}}\frac{\partial \widehat{T}}{\partial \widehat{y}}+\frac{D_T}{{\widehat{T}}_{\infty }}{\left(\frac{\partial \widehat{T}}{\partial \widehat{y}}\right)}^2\right\}-\frac{1}{{\left(\rho C_P\right)}_P}\frac{\partial q_r}{\partial \widehat{y}}$$

(3)

Boundary layer form of mass concentration equation

$$\widehat{u}\frac{\partial \widehat{C}}{\partial \widehat{x}}+\widehat{v}\frac{\partial \widehat{C}}{\partial \widehat{y}}=D_B\frac{{\partial }^2\widehat{C}}{{\partial }^2\widehat{y}}+\frac{D_T}{T_{\infty }}\frac{{\partial }^2\widehat{T}}{{\partial }^2\widehat{y}}$$

(4)

Boundary conditions

$$\begin{gathered} \widehat{u}=0=\widehat{v}, \widehat{T}={\widehat{T}}_w, \widehat{C}={\widehat{C}}_{w}at\widehat{y}=0 \\ \widehat{u}\to 0=\widehat{v}, \widehat{T}\to {\widehat{T}}_{\infty }, \widehat{C}\to {\widehat{C}}_{\infty }as\widehat{y}\to \infty \end{gathered}$$

(5)

The equation of continuity (Equation (1)) is also called the law of conservation of the mass; it describes that the amount of mass entering per unit area is equal to the amount of mass leaving per unit area. Alternatively, we can say that mass remains conserved. The momentum equation (Equation (2)) is also called the approximate form of the boundary layer Navier–Stokes equation. It has derived from Newton’s second law of motion. This equation describes the motion of the fluid. Equation (3) corresponds to the energy equation—it is also called the law of conservation of energy, and it has derived from the first law thermodynamics. This equation states that the entire energy of the system remains conserved. Equation (4) is called mass concentration and the law conservation of concentration. Equation (5) represents the imposed boundary conditions for the current molded process. Here, $\widehat{u}$ and $\widehat{v}$ are the velocity components along $\widehat{x}$ and $\widehat{y}$ respectively. The variables $r\left(\widehat{x}\right)=asin\widehat{x}$ are the radial distance from the symmetric axis to surface of the sphere. The designations $\nu =\frac{\mu }{\rho },g,{\beta }_T,$ and ${\beta }_C$ represent the kinematic viscosity with $\mu$ as fluid dynamic viscosity, gravitation acceleration, expansion coefficient due to temperature, and expansion coefficient due to concentration, respectively. Notations ${\alpha }_m=k/{\left(\rho C_P\right)}_P,D_B,D_T,\rho$ and $C_P$ denote thermal diffusivity with $k$ as the thermal conductivity of the fluid, Brownian movement coefficient, thermophoresis coefficient, density of the fluid, and specific heat at constant pressure, respectively. Symbols ${\left(\rho C_P\right)}_P,{\left(\rho C_P\right)}_f,$ and $\tau =\frac{{\left(\rho C_P\right)}_P}{{\left(\rho C_P\right)}_f}$ denote heat capacity of the nanofluid, heat capacity of the base fluid, and the ratio between the heat capacities, respectively. The notation $(T,T_{\infty })$ represents liquid and ambient temperatures, and $(C,C_{\infty })$ represents liquid and ambient concentrations.

Using the Rosseland approximation [37,38,39] for radiation in an optically thick fluid, the radiative heat flux is simplified to:

$$q_r=-\frac{4\sigma }{3k^{*}}\frac{\partial T^4}{\partial y}$$

(6)

where $\sigma$ and $k^{*}$ designates the Stefan–Boltzmann constant and mean absorption coefficient, respectively. We assume that the temperature differences within the flow, such as the term $T^4$, may be expressed as a linear function of temperature. Hence, expanding $T^4$ in a Taylor series about $T_{\infty }$ and neglecting higher-order terms, we obtain

$$T^4=4T_{\infty }T-3{T_{\infty }}^4$$

(7)

In view of Equations (6) and (7), Equation (3) takes the following form:

$$\widehat{u}\frac{\partial \widehat{T}}{\partial \widehat{x}}+\widehat{v}\frac{\partial \widehat{T}}{\partial \widehat{y}}={\alpha }_m\left(1+\frac{4}{3N}\right)\frac{{\partial }^2\widehat{T}}{{\partial }^2\widehat{y}}+\tau \left\{D_B\frac{\partial \widehat{C}}{\partial \widehat{y}}\frac{\partial \widehat{T}}{\partial \widehat{y}}+\frac{D_T}{{\widehat{T}}_{\infty }}{\left(\frac{\partial \widehat{T}}{\partial \widehat{y}}\right)}^2\right\}$$

(8)

where $N=\frac{\rho C_{P}k{k}^{*}}{4T_{\infty }^3\sigma }$ is the radiation parameter.

3. Solution Process

Now, in this section, we elaborate the entire solution procedure for the flow model presented in Equations (1)–(4) with boundary conditions given in Equation(5).

3.1. Dimensionless Variables

In the current study, the finite difference method is applied to the set of partial differential equations to determine the numerical solutions. Before the application of finite difference method, the model for the natural convection along the surface of the sphere and in the plume-region is made dimensionless (unit-free) by using the suitable dimensionless variables defined in the following Equation (9) used in Refs. [8,9]

$$\begin{gathered} x=\frac{\widehat{x}}{a},y={Gr}^{\frac{1}{4}}\frac{\widehat{y}}{y},u={G_r}^{-\frac{1}{2}}\frac{a}{\nu }\widehat{u},v={Gr}^{-\frac{1}{4}}\frac{a}{\nu }\widehat{v},r=\frac{\widehat{r}}{a},\theta =\frac{{T-T}_{\infty }}{T_w-T_{\infty }},\varphi =\frac{{C-C}_{\infty }}{C_w-C_{\infty }}, \\ G{r}_T=\frac{g\beta ∆T{a}^3}{{\nu }^2},G{r}_C=\frac{g{\beta }_c∆C{a}^3}{{\nu }^2},Gr=G{r}_T+G{r}_C. \end{gathered}$$

(9)

where ${Gr}_T=\frac{g\beta ∆\widehat{T}{a}^3}{{\nu }^2}$ is the Grashoff number,$G{r}_C=\frac{g{\beta }_c∆T{a}^3}{{\nu }^2}$ is the modified Grashof number, $\theta$ is the dimensionless temperature function, $\varphi$ is the dimensionless mass concentration, $x$ is the dimensionless horizontal axis and $y$ is the dimensionless vertical axis, $a$ is the radius of a sphere, $u$ is the non-dimensional horizontal velocity component, and $v$ is the non-dimensional velocity component. By using the dimensionless variables in Equation (9) into Equations (1)–(4) along with boundary condition (5), we have the following dimensionless governing equations:

Dimensionless form of conservation of mass equation

$$\frac{\partial \left(usinx\right)}{\partial x}+\frac{\partial \left(vsinx\right)}{\partial y}=0$$

(10)

Dimensionless form of momentum equation

$$u\frac{\partial u}{\partial x}+v\frac{\partial u}{\partial y}=\frac{{\partial }^{2}u}{{\partial }^{2}y}+\theta Sinx+\varphi Sinx$$

(11)

Dimensionless form of energy equation

$$u\frac{\partial \theta }{\partial x}+v\frac{\partial \theta }{\partial y}=\frac{1}{P_r}\left[1+\frac{4}{3N}\right]\frac{{\partial }^2\theta }{{\partial }^{2}y}+N_b\frac{\partial \phi }{\partial y}\frac{\partial \theta }{\partial y}+N_t{\left(\frac{\partial \theta }{\partial y}\right)}^2$$

(12)

Dimensionless form of momentum equation

$$u\frac{\partial \varphi }{\partial x}+v\frac{\partial \varphi }{\partial y}=\frac{1}{Sc}\left(\frac{{\partial }^2\phi }{{\partial }^{2}y}+\frac{N_t}{N_b}\frac{{\partial }^2\theta }{{\partial }^{2}y}\right)$$

(13)

Dimensionless form subjected boundary conditions,

$$\begin{gathered} u=0=v, \varphi =1,\hspace{1em}\theta =1\hspace{1em}aty=0 \\ u\to 0,\hspace{1em}\hspace{1em}\varphi \to 0,\hspace{1em}\theta \to 0\hspace{1em}asy\to \infty \end{gathered}$$

(14)

Here, the physical parameter $Pr=\frac{\nu }{{\alpha }_m}$ is the Prandtl number, $Sc=\frac{\nu }{D_B}$ is the Schmidt number, $N_b=\frac{{\left(\rho c\right)}_p{D}_B\left(C_w-C_{\mathrm{\infty }}\right)}{{\left(\rho c\right)}_f\nu }$ is the Brownian motion parameter,$N_t=\frac{{\left(\rho c\right)}_p{D}_T\left(T_w-T_{\mathrm{\infty }}\right)}{{\left(\rho c\right)}_f{T}_{\mathrm{\infty }}\nu }$ is the thermophoresis parameter, and $N=\frac{\rho C_{P}k{k}^{\mathrm{*}}}{4T_{\mathrm{\infty }}^3\sigma }$ is the radiation parameter.

Moreover, the skin friction coefficient, Nusselt number, and Sherwood number at the surface can be found as follows:

$$C_f=\frac{{a\tau }_w}{\rho U^2},Nu=\frac{a{q}_w}{k(T_w-T_{\infty })},Sh=\frac{a{J}_w}{D_B(C_w-C_{\infty })}$$

(15)

where:

$${\tau }_w={\mu \left(\frac{\partial \stackrel{-}{u}}{\partial \stackrel{-}{y}}\right)}_{\stackrel{-}{y}=0},q_w={-k\left(\frac{\partial T}{\partial \stackrel{-}{y}}\right)}_{\stackrel{-}{y}=0},J_w={-D_B\left(\frac{\partial C}{\partial \stackrel{-}{y}}\right)}_{\stackrel{-}{y}=0}$$

(16)

By using Equation (9) into Equations (14) and (15), we have the following expressions:

$$C_{f}G{r}^{1/4}={\left(\frac{\partial u}{\partial y}\right)}_{y=0}$$

(17)

$$NuG{r}^{-1/4}={-\left(\frac{\partial \theta }{\partial y}\right)}_{y=0}$$

(18)

$$ShG{r}^{-1/4}={-\left(\frac{\partial \varphi }{\partial y}\right)}_{y=0}$$

(19)

3.2. Primitive Variable Formulation

To reduce the labor on the numerical algorithm for the finite difference method on FORTRAN Lahy-95, the dimensionless model presented in Equations (10)–(13) along with boundary conditions (14) is further transformed to a primitive form by using the following appropriate primitive variables formulation used by [8,9]:

$$u=x^{\frac{1}{2}}U\left(X,Y\right),v=x^{-\frac{1}{4}}V\left(X,Y\right),Y=x^{-\frac{1}{4}}yx=X,\theta =\stackrel{-}{\theta }\left(X,Y\right),\hspace{1em}\varphi =\stackrel{-}{\varphi }\left(X,Y\right),r\left(x\right)=asinX$$

(20)

After using Equation (20) into Equations (10)–(14), the following equations are attained:

$$X^{1/2}UcosX+\left\{X\frac{\partial U}{\partial X}-\frac{1}{4}Y\frac{\partial U}{\partial Y}+\frac{1}{2}U+\frac{\partial V}{\partial Y}\right\}sinX=0$$

(21)

$$XU\frac{\partial U}{\partial X}+\frac{1}{2}{U}^2+\left(V-\frac{1}{4}YU\right)\frac{\partial U}{\partial Y}=\frac{{\partial }^{2}U}{{\partial Y}^2}+\stackrel{-}{\theta }sinX+\stackrel{-}{\phi }sinX$$

(22)

$$XU\frac{\partial \stackrel{-}{\theta }}{\partial X}+\left(V-\frac{1}{4}YU\right)\frac{\partial \stackrel{-}{\theta }}{\partial Y}=\frac{1}{Pr}\left[1+\frac{4}{3N}\right]\frac{{\partial }^2\stackrel{-}{\theta }}{{\partial Y}^2}+N_b\frac{\partial \stackrel{-}{\varphi }}{\partial Y}\frac{\partial \stackrel{-}{\theta }}{\partial Y}$$

(23)

$$XU\frac{\partial \stackrel{-}{\varphi }}{\partial X}+\left(V-\frac{1}{4}YU\right)\frac{\partial \stackrel{-}{\varphi }}{\partial Y}=\frac{1}{Sc}\left(\frac{{\partial }^2\stackrel{-}{\varphi }}{{\partial }^{2}Y}+\frac{N_t}{N_b}\frac{{\partial }^2\stackrel{-}{\theta }}{{\partial }^{2}Y}\right)$$

(24)

The subjected conditions are:

$$\begin{gathered} U=0=V,\stackrel{-}{\varphi }=1,\stackrel{-}{\theta }=1\hspace{1em}atY=0 \\ U\to 0,\stackrel{-}{\varphi }\to 0,\stackrel{-}{\theta }\to 0\hspace{1em}asY\to \infty . \end{gathered}$$

(25)

The transformed model in primitive form is presented in Equations (21)–(24) along with the transformed boundary conditions given in Equation (25). The transformed boundary conditions are so smooth and have similar terms which make coding of the algorithm easy.

3.3. Computational Technique

A finite difference method is taken into use for numerical solutions of altered flow equations, as shown in Equations (21)–(24), with the boundary conditions shown in Equation (25). The $X$-axis is used to apply the backward difference, while the $Y$-axis is used to apply the central difference. After the flow equations are discretized, the unknown discretized variables ${(U}_{i,j},V_{i,j},{\theta }_{i,j},{\mathrm{\varphi }}_{\mathrm{I},\mathrm{j}})$ are identified. The finite difference approach is used to arrive at the solution of Equations (21)–(24) with problem conditions specified in Equation (25).The backward and center differences are applied in $x$ and $y$ directions, respectively. The following is a description of the solution process:

$$\frac{\partial U}{\partial X}=\frac{U_{(i,j)}-U_{(i,j-1)}}{∆X}$$

(26)

$$\frac{\partial U}{\partial Y}=\frac{U_{(i+1,j)}-U_{(i,-1,j)}}{2∆Y}$$

(27)

$$\frac{{\partial }^{2}U}{\partial Y^2}=\frac{U_{(i+1,j)}-2U_{(i,j)}+U_{(i,-1,j)}}{∆Y^2}$$

(28)

The following system of algebraic equations is created by incorporating Equations (26)–(28), Equations (21)–(24), and the boundary conditions specified in Equation (25).

The discretized form of continuity equation,

$$V_{i,j}=\left\{V_{i-1,j}-X_i\frac{∆Y}{∆X}\left(U_{i,j}-U_{i,j-1}\right)+\frac{1}{8}{Y}_j\left(U_{i+1,j}-U_{i-1,j}\right)-\frac{1}{2}{∆YU}_{i,j}\right\}-\frac{cos{X}_i}{sin{X}_i}{{X_i}^{\frac{1}{2}}U}_{i,j}$$

(29)

The discretized form of momentum equation,

$$\begin{gathered} \left(\frac{1}{2}∆Y\left(V_{i,j}-\frac{1}{4}{Y}_j{U}_{i,j}\right)+1\right)U_{i-1,j}+\left(-{∆Y}^2\left(\frac{1}{2}+\frac{X_i}{∆X}\right)U_{i,j}-2\right)U_{i,j}+\left(-\frac{1}{2}∆Y\left(V_{i,j}-\frac{1}{4}{Y}_j{U}_{i,j}\right)+1\right){\mathit{U}}_{\mathit{i}+1,\mathit{j}} \\ ={-X}_i{U}_{i,j}{U}_{i,j-1}\frac{{∆Y}^2}{∆X}-{∆Y}^2\left({\stackrel{-}{\theta }}_{i,j}sin{X}_i+{\stackrel{-}{\varphi }}_{i,j}sin{X}_i\right) \end{gathered}$$

(30)

The discretized form of energy equation,

$$\begin{gathered} \left(\frac{1}{2}∆Y\left(V_{i,j}-\frac{1}{4}{Y}_j{U}_{i,j}\right)+\frac{1}{P_r}\left[1+\frac{4}{3N}\right]-N_b\frac{1}{4}\left({\stackrel{-}{\varphi }}_{i+1,j}-{\stackrel{-}{\varphi }}_{i-1,j}\right)+N_t\frac{1}{4}{\stackrel{-}{\theta }}_{i-1,j}\right){\stackrel{-}{\theta }}_{i-,j}+ \\ \left(-X_i{U}_{i,j}\frac{{∆Y}^2}{∆X}-\frac{2}{P_r}\left[1+\frac{4}{3N}\right]\right){\stackrel{-}{\theta }}_{i,j}+\left(-\frac{1}{2}∆Y\left(V_{i,j}-\frac{1}{4}{Y}_j{U}_{i,j}\right)+\frac{1}{P_r}\left[1+\frac{4}{3N}\right] \\ +\frac{1}{4}\left({\stackrel{-}{\phi }}_{i+1,j}-{\stackrel{-}{\phi }}_{i-1,j}\right)+N_t\frac{1}{4}{\stackrel{-}{\theta }}_{i+1,j}\right){\stackrel{-}{\theta }}_{i+,j}={-X}_i{U}_{i,j}\frac{{∆Y}^2}{∆X}{\stackrel{-}{\theta }}_{i,j-1}+N_t\frac{1}{2}{\stackrel{-}{\theta }}_{i+1,j}{\stackrel{-}{\theta }}_{i-1,j} \end{gathered}$$

(31)

The discretized form of mass concentration equation,

$$\begin{gathered} \left(\frac{1}{2}∆Y\left(V_{i,j}-\frac{1}{4}{Y}_j{U}_{i,j}\right)+\frac{1}{Sc}\right){\stackrel{-}{\varphi }}_{i-1,j}+\left({-X}_i{U}_{i,j}\frac{{∆Y}^2}{∆X}-\frac{2}{Sc}\right){\stackrel{-}{\varphi }}_{\mathit{i},\mathit{j}}+\left(-\frac{1}{2}∆Y\left(V_{i,j}-\frac{1}{4}{Y}_j{U}_{i,j}\right)+\frac{1}{Sc}\right){\stackrel{-}{\varphi }}_{\mathit{i}+1,\mathit{j}} \\ =-X_i{U}_{i,j}\frac{{∆Y}^2}{∆X}{\stackrel{-}{\theta }}_{i,j-1}+\frac{1}{Sc}\frac{N_t}{N_b}\left({\stackrel{-}{\theta }}_{i+1,j}-2{\stackrel{-}{\theta }}_{i,j}+{\stackrel{-}{\theta }}_{i-1,j}\right) \end{gathered}$$

(32)

The subjected, discretized boundary conditions,

$$\begin{gathered} U_{i,j}=0=V_{i,j},\hspace{1em}{\stackrel{-}{\theta }}_{i,j}=1={\stackrel{-}{\varphi }}_{i,j}\hspace{1em}aty=0 \\ {U}_{i,j}\to 0,\hspace{1em}{V}_{i,j}\to 0,{\stackrel{-}{\theta }}_{i,j}\to 0,{\stackrel{-}{\varphi }}_{i,j}\to 0\hspace{1em}asy\to \infty \end{gathered}$$

(33)

The approximate solutions determined with the finite difference method are discussed in detail in the forthcoming section. The convergence and mesh sensitivity analysis criterion are presented as follows to achieve accurate numerical solutions for the variables $U,V$, $\mathrm{\theta }$ and $\varphi$, respectively,

$$\max \left|U_{ij}\right|+\max \left|V_{ij}\right|+\max \left|{\mathsf{\theta }}_{ij}\right|+\max \left|{\mathsf{\varphi }}_{ij}\right|\le ϵ$$

where $ϵ={10}^{-5}$. The computation is started at $X=0$ and then marches downstream implicitly. Here, we have taken the step size $\mathrm{\Delta }X=0.05$ and $\mathrm{\Delta }Y=0.02$, respectively.

4. Flow Analysis for Plume Region-III

The plume region is established just above the sphere due to eruption fluid from the boundary layer around the sphere. The thermal radiation effects along with thermophoresis and Brownian motion effects are considered in the following model. The base fluid for the making of nanofluids is water, and due to radiative heat flux, it is optically dens gray. The surface of the sphere is maintained by the temperature ${\widehat{T}}_w$ and its temperature is considered greater than the ambient temperature ${\widehat{T}}_{\infty }$. The mass concentration on surface of the sphere is ${\widehat{C}}_w$ and this concentration is larger than the concentration in the free-stream region ${\widehat{C}}_{\infty }$. The axes $\widehat{x}$ and $\widehat{z}$ are taken for the flow direction and normal to this flow, respectively; this coordinated system is highlighted in Figure 1. The velocity components corresponding to the flow direction and normal to this flow are designated by $\widehat{u}$ and $\widehat{W}$, respectively. Three regions make up the flow’s schematic diagram: region-I is the spherical boundary layer, region-II is the eruption region, and region-III is the plume in Figure 1. Following [1,2,3,4,8,9], the fundamental flow equations are shown below:

The law of conservation of mass (equation of continuity),

$$\frac{\partial \left(\widehat{z}\widehat{u}\right)}{\partial \widehat{x}}+\frac{\partial \left(\widehat{z}\widehat{w}\right)}{\partial \widehat{z}}=0$$

(34)

Boundary layer form of momentum equation

$$\widehat{u}\frac{\partial \widehat{u}}{\partial \widehat{x}}+\widehat{w}\frac{\partial \widehat{u}}{\partial \widehat{z}}=\nu \frac{1}{\widehat{z}}\frac{\partial }{\partial \widehat{z}}\left(\widehat{z}\frac{\partial \widehat{u}}{\partial \widehat{z}}\right)-g\beta \left(\widehat{T}-{\widehat{T}}_{\infty }\right)-g{\beta }_c\left(\widehat{C}-{\widehat{C}}_{\infty }\right)$$

(35)

Boundary layer form of energy equation

$$\widehat{u}\frac{\partial \widehat{T}}{\partial \widehat{x}}+\widehat{w}\frac{\partial \widehat{T}}{\partial \widehat{z}}=\alpha \frac{1}{\widehat{z}}\frac{\partial }{\partial \stackrel{-}{z}}\left(\stackrel{-}{z}\frac{\partial \widehat{T}}{\partial \widehat{z}}\right)+\tau \left\{D_B\frac{\partial \stackrel{-}{C}}{\partial \widehat{z}}\frac{\partial \widehat{T}}{\partial \widehat{z}}+\frac{D_T}{T_{\infty }}{\left(\frac{\partial \widehat{T}}{\partial \widehat{z}}\right)}^2\right\}-\frac{1}{\rho c_p}\frac{\partial q_r}{\partial \widehat{z}}$$

(36)

Boundary layer form of mass concentration equation

$$\widehat{u}\frac{\partial \widehat{C}}{\partial \widehat{x}}+\widehat{w}\frac{\partial \widehat{C}}{\partial \widehat{z}}=D_B\frac{1}{\widehat{z}}\frac{\partial }{\partial \widehat{z}}\left(\widehat{z}\frac{\partial \widehat{C}}{\partial \widehat{z}}\right)+\frac{D_T}{{\widehat{T}}_{\infty }}\frac{1}{\widehat{z}}\frac{\partial }{\partial \widehat{z}}\left(\widehat{z}\frac{\partial \widehat{T}}{\partial \widehat{z}}\right)$$

(37)

Molded conditions

$$\begin{gathered} \widehat{u}=\widehat{w}=0,\widehat{T}={\widehat{T}}_w,\widehat{C}={\widehat{C}}_{w}at\widehat{z}=0 \\ \widehat{u}\to 0,\widehat{T}\to {\widehat{T}}_{\infty },\hspace{1em}\hspace{1em}\widehat{C}\to {\widehat{C}}_{\infty }as\widehat{z}\to \infty \end{gathered}$$

(38)

The equation of continuity (Equation (34)) is also called the law of conservation of mass for plume region-III—it describes that the amount of mass entering per unit area is equal to the amount of mass leaving per unit area. Alternatively, we can say that mass remains conserved. The momentum equation (Equation (35)) is also called the approximate form of the boundary layer Navier–Stokes equation for plume region-III. It has derived from Newton’s second law of motion. This equation describes the motion of the fluid. Equation (36) is called the energy equation and the law of conservation of energy; it has derived from the first law thermodynamics for plume region-III. This equation states that the entire energy of the system remains conserved. Equation (37) is called mass concentration and is called the law conservation of concentration for plume region-III. Equation (38) represents the imposed boundary conditions for the current molded process for plume region-III. Here, $\widehat{u}$ and $\widehat{v}$ are the velocity components along $\widehat{x}$ and $\widehat{y}$, respectively. The variable $r\left(\widehat{x}\right)=asin\widehat{x}$ is the radial distance from the symmetric axis to the surface of the sphere. The designations $\nu =\frac{\mu }{\rho },g,{\beta }_T,$ and ${\beta }_C$ represent the kinematic viscosity with $\mu$ as fluid dynamic viscosity, gravitation acceleration, expansion coefficient due to temperature, and expansion coefficient due to concentration, respectively. Notations ${\alpha }_m=k/{\left(\rho C_P\right)}_P,D_B,D_T,\rho$ and $C_P$ denote thermal diffusivity with $k$ as the thermal conductivity of the fluid, Brownian movement coefficient, thermophoresis coefficient, density of the fluid, and specific heat at a constant pressure, respectively. Symbols ${\left(\rho C_P\right)}_P,{\left(\rho C_P\right)}_f,$ and $\tau =\frac{{\left(\rho C_P\right)}_P}{{\left(\rho C_P\right)}_f}$ denote the heat capacity of the nanofluid, heat capacity of the base fluid, and the ratio between the heat capacities, respectively. The notation $(T,T_{\infty })$ represents liquid and ambient temperatures, and $\left(C,C_{\infty }\right)$ represents liquid and ambient concentrations.

Using the Rosseland approximation for radiation in an optically thick fluid, [37,38,39] the radiative heat flux is simplified to:

$$q_r=-\frac{4\sigma }{3k^{*}}\frac{\partial T^4}{\partial y}$$

(39)

where $\sigma$ and $k^{*}$ designate the Stefan–Boltzmann constant and the mean absorption coefficient, respectively. We assume that the temperature differences within the flow, such as the term $T^4$, may be expressed as a linear function of temperature. Hence, by expanding $T^4$ in a Taylor series about $T_{\infty }$ and neglecting higher-order terms, we obtain

$$T^4=4T_{\infty }T-3{T_{\infty }}^4$$

(40)

In view of Equation (6) and (7), Equation (3) takes the following form:

$$\widehat{u}\frac{\partial \widehat{T}}{\partial \widehat{x}}+\widehat{v}\frac{\partial \widehat{T}}{\partial \widehat{y}}={\alpha }_m\left(1+\frac{4}{3N}\right)\frac{{\partial }^2\widehat{T}}{{\partial }^2\widehat{y}}+\tau \left\{D_B\frac{\partial \widehat{C}}{\partial \widehat{y}}\frac{\partial \widehat{T}}{\partial \widehat{y}}+\frac{D_T}{{\widehat{T}}_{\infty }}{\left(\frac{\partial \widehat{T}}{\partial \widehat{y}}\right)}^2\right\}$$

(41)

where $N=\frac{\rho C_{P}k{k}^{*}}{4T_{\infty }^3\sigma }$ is the is the radiation parameter.

5. Solution Process

Now, in this section, we elaborate the entire solution procedure for the flow model presented in Equations (34)–(37) with the boundary conditions given in Equation(38).

5.1. Dimensionless Variables

In the current study, the finite difference method is applied on the system of partial differential equations to figure out the numerical solutions. Before the application of finite difference method, the model for the natural convection along the surface of the sphere and in the plume-region is made dimensionless (unit-free) by using the suitable dimensionless variables defined in the following Equation (39) used in Refs. [8,9]

$$\begin{gathered} x=\frac{\widehat{x}}{a},z={Gr}^{\frac{1}{4}}\frac{\widehat{z}}{z},u={Gr}^{-\frac{1}{2}}\frac{a}{\nu }\widehat{u},w={Gr}^{-\frac{1}{4}}\frac{a}{\nu }\widehat{w},r=\frac{\widehat{r}}{a},\theta =\frac{{T-T}_{\infty }}{T_w-T_{\infty }},\varphi =\frac{{C-C}_{\infty }}{C_w-C_{\infty }}, \\ G{r}_T=\frac{g\beta ∆T{a}^3}{{\nu }^2},G{r}_C=\frac{g{\beta }_c∆C{a}^3}{{\nu }^2},Gr=G{r}_T+G{r}_C. \end{gathered}$$

(42)

where $G_T=\frac{g\beta ∆\widehat{T}{a}^3}{{\nu }^2}$ is the Grashoff number,$G{r}_C=\frac{g{\beta }_c∆C{a}^3}{{\nu }^2}$ is the changed Grashof number, $\theta$ is the dimensionless temperature function, $\varphi$ is the dimensionless mass concentration, $x$ is the dimensionless horizontal axis and $z$ is the dimensionless vertical axis, $a$ is the radius of a sphere, $u$ is the non-dimensional horizontal velocity component, and $w$ is the non-dimensional velocity component. By using the dimensionless variables in Equation (42) into Equations (34)–(37) along with boundary condition (38), we have the following dimensionless governing equations:

Non-dimensional law of conservation of mass (equation of continuity),

$$\frac{\partial \left(ru\right)}{\partial x}+\frac{\partial \left(rw\right)}{\partial z}=0,$$

(43)

Non-dimensional boundary layer form of the momentum equation

$$u\frac{\partial u}{\partial x}+w\frac{\partial u}{\partial z}=\frac{1}{z}\frac{\partial }{\partial z}\left(z\frac{\partial u}{\partial z}\right)-\theta -\varphi ,$$

(44)

Non-dimensional boundary layer form of the energy equation

$$u\frac{\partial \theta }{\partial x}+w\frac{\partial \theta }{\partial z}=\frac{1}{Pr}\frac{1}{z}\frac{\partial }{\partial z}\left(z\frac{\partial \theta }{\partial z}\right)+N_b\frac{\partial \phi }{\partial z}\frac{\partial \theta }{\partial z}+N_t{\left(\frac{\partial \theta }{\partial z}\right)}^2+\frac{4}{3NPr}\frac{{\partial }^2\theta }{{\partial }^{2}z},$$

(45)

Non-dimensional boundary layer form of the mass concentration equation

$$u\frac{\partial \varphi }{\partial x}+w\frac{\partial \varphi }{\partial z}=\frac{1}{Sc}\left(\frac{1}{z}\frac{\partial }{\partial z}\left(z\frac{\partial \phi }{\partial z}\right)+\frac{N_t}{N_b}\frac{1}{z}\frac{\partial }{\partial z}\left(z\frac{\partial \theta }{\partial z}\right)\right),$$

(46)

Non-dimensional appropriate boundary conditions

$$\begin{gathered} u=0,w=0 \theta =1,\hspace{1em}\hspace{1em}\varphi =1\hspace{1em}atz=0 \\ u\to 0,\hspace{1em}\hspace{1em}\varphi \to 0,\theta \to 0\hspace{1em}\hspace{1em}asz\to \infty \end{gathered}$$

(47)

Here, the physical parameters $Pr=\frac{\nu }{{\alpha }_m}$ is the Prandtl number, $Sc=\frac{\nu }{D_B}$ is the Schmidt number, $N_b=\frac{{\left(\rho c\right)}_p{D}_B\left(C_w-C_{\mathrm{\infty }}\right)}{{\left(\rho c\right)}_f\nu }$ is the Brownian motion parameter, $N_t=\frac{{\left(\rho c\right)}_p{D}_T\left(T_w-T_{\mathrm{\infty }}\right)}{{\left(\rho c\right)}_f{T}_{\mathrm{\infty }}\nu }$ is the thermophoresis parameter, and $N=\frac{\rho C_{P}k{k}^{\mathrm{*}}}{4T_{\mathrm{\infty }}^3\sigma }$ is the radiation parameter.

Moreover, the skin friction coefficient, Nusselt number, and Sherwood number at the surface can be found as follows:

$$C_f=\frac{{a\tau }_w}{\rho U^2},Nu=\frac{a{q}_w}{k(T_w-T_{\infty })},Sh=\frac{a{J}_w}{D_B(C_w-C_{\infty })}$$

(48)

where:

$${\tau }_w={\mu \left(\frac{\partial u}{\partial z}\right)}_{\stackrel{-}{y}=0},q_w={-k\left(\frac{\partial T}{\partial z}\right)}_{\stackrel{-}{y}=0},J_w={-D_B\left(\frac{\partial C}{\partial z}\right)}_{\stackrel{-}{y}=0}$$

(49)

By using Equation (42) into Equations (48)–(49), we have the following expressions:

$$C_{f}G{r}^{1/4}={\left(\frac{\partial u}{\partial z}\right)}_{y=0}$$

(50)

$$NuG{r}^{-1/4}={-\left(\frac{\partial \theta }{\partial z}\right)}_{y=0}$$

(51)

$$ShG{r}^{-1/4}={-\left(\frac{\partial \varphi }{\partial z}\right)}_{y=0}$$

(52)

5.2. Primitive Variable Formulation

To simplify the numerical algorithm for the finite difference method on FORTRAN Lahy-95, the dimensionless model presented in Equations (43)–(46) along with boundary conditions (47) is further transformed to a primitive form by using the following appropriate primitive variable formulation used by [8,9]:

$$u=x^{\frac{1}{2}}U\left(X,Z\right),W=x^{-\frac{1}{4}}W\left(X,Z\right),x=X,Z=x^{-1/4}Z,\hspace{1em}\hspace{1em}\theta =\stackrel{-}{\theta }\left(X,Z\right),\hspace{1em}\varphi =\stackrel{-}{\varphi }\left(X,Z\right)$$

(53)

By using the variables given in Equation (53) into Equations (43)–(46) along with the boundary conditions given in Equation(47), we have the following expressions:

$$Z\frac{\partial U}{\partial X}-\frac{Z^2}{4X}\frac{\partial U}{\partial Z}+\frac{3}{4}ZU+W+Z\frac{\partial W}{\partial Z}=0$$

(54)

$$XU\frac{\partial U}{\partial X}+\frac{1}{2}{U}^2+\left(W-\frac{1}{4}ZU\right)\frac{\partial U}{\partial Z}=\frac{1}{Z}\frac{\partial }{\partial Z}\left(Z\frac{\partial U}{\partial Z}\right)-\stackrel{-}{\theta }-\stackrel{-}{\varphi }$$

(55)

$$XU\frac{\partial \stackrel{-}{\theta }}{\partial X}+\left(W-\frac{1}{4}ZU\right)\frac{\partial \stackrel{-}{\theta }}{\partial Z}=\frac{1}{P_r}\frac{1}{Z}\frac{\partial }{\partial Z}\left(Z\frac{\partial \stackrel{-}{\theta }}{\partial Z}\right)+N_b\frac{\partial \stackrel{-}{\varphi }}{\partial Z}\frac{\partial \stackrel{-}{\theta }}{\partial Z}+N_t{\left(\frac{\partial \stackrel{-}{\theta }}{\partial Z}\right)}^2+\frac{4}{3N{P}_r}\frac{{\partial }^2\stackrel{-}{\theta }}{\partial Z^2}$$

(56)

$$XU\frac{\partial \stackrel{-}{\varphi }}{\partial X}+\left(W-\frac{1}{4}ZU\right)\frac{\partial \stackrel{-}{\varphi }}{\partial Z}=\frac{1}{Sc}\left(\frac{1}{Z}\frac{\partial }{\partial Z}\left(Z\frac{\partial \stackrel{-}{\varphi }}{\partial Z}\right)+\frac{N_t}{N_b}\frac{1}{Z}\frac{\partial }{\partial Z}\left(Z\frac{\partial \stackrel{-}{\theta }}{\partial Z}\right)\right)$$

(57)

Boundary conditions

$$\begin{gathered} W=U=0,\hspace{1em}\stackrel{-}{\varphi }=1,\hspace{1em}\stackrel{-}{\theta }=1\hspace{1em}atZ=0 \\ U\to 0,\hspace{1em}\hspace{1em}\stackrel{-}{\varphi }\to 0, \hspace{1em}\stackrel{-}{\theta }\to 0\hspace{1em} asZ\to \infty \end{gathered}$$

(58)

The transformed model in primitive form is presented in Equations (54)–(57) along with the transformed boundary conditions given in Equation (58); the transformed boundary conditions are so smooth and have similar terms which make coding of the algorithm easy.

5.3. Computational Technique

A finite difference method is taken into use for the numerical solutions of altered flow equations, as shown in Equations (21)–(24), with the boundary conditions shown in Equation (25). The $X$-axis is used to apply the backward difference, while the $Y$-axis is used to apply the central difference. After the flow equations are discretized, the unknown discretized variables ${(U}_{i,j},W_{i,j},{\theta }_{i,j},{\mathrm{\varphi }}_{\mathrm{i},\mathrm{j}})$ are identified. The finite difference approach is used to arrive at the solution of Equations (54)–(57) with the problem conditions specified in Equation (58).The backward and center differences are applied in $x$ and $y$ directions, respectively. The following is a description of the solution process:

$$\frac{\partial U}{\partial X}=\frac{U_{(i,j)}-U_{(i,j-1)}}{∆X}$$

(59)

$$\frac{\partial U}{\partial Z}=\frac{U_{(i+1,j)}-U_{(i,-1,j)}}{2∆Z}$$

(60)

$$\frac{{\partial }^{2}U}{\partial Z^2}=\frac{U_{(i+1,j)}-2U_{(i,j)}+U_{(i,-1,j)}}{∆Z^2}$$

(61)

The following system of algebraic equations is created by incorporating Equations (59)–(61),Equations (54)–(57), and the boundary conditions specified in Equation (58).

$$W_{i,j}=\frac{1}{\left(∆Z+Z_j\right)}\left\{Z_j{W}_{i-1,j}-Z_j\left(U_{i,j}-U_{i,j-1}\right)\frac{∆Z}{∆X}+\frac{{Z_j}^2}{8X_i}\left(U_{i+1,j}-U_{i-1,j}\right)-\frac{3}{4}{∆ZZ}_j{U}_{i,j}\right\}$$

(62)

$$\begin{gathered} \left(\left[\frac{1}{2}∆Z\left(W_{i,j}-\frac{1}{4}{Z}_j{U}_{i,j}\right)-\frac{∆Z}{{2Z}_j}\right]+1\right){\mathit{U}}_{\mathit{i}-1,\mathit{j}}+\left(-\left(X_i\frac{{∆Z}^2}{∆X}+\frac{1}{2}{∆Z}^2\right)U_{i,j}-2\right){\mathit{U}}_{\mathit{i},\mathit{j}} \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+\left(-\left[\frac{1}{2}∆Z\left(W_{i,j}-\frac{1}{4}{Z}_j{U}_{i,j}\right)-\frac{∆Z}{{2Z}_j}\right]+1\right){\mathit{U}}_{\mathit{i}-1,\mathit{j}} \\ =-{∆Z}^2\frac{X_i}{∆X}{U}_{i,j}{U}_{i,j-1}+{\stackrel{-}{\theta }}_{i,j}+{\stackrel{-}{\phi }}_{i,j} \end{gathered}$$

(63)

$$\begin{gathered} \left(\left[\frac{1}{2}∆Z\left(W_{i,j}-\frac{1}{4}{Z}_j{U}_{i,j}\right)-\frac{1}{P_r}\frac{∆Z}{2Z_j}-N_b\frac{1}{4}\left({\stackrel{-}{\phi }}_{i+1,j}-{\stackrel{-}{\phi }}_{i-1,j}\right)+N_t\frac{1}{4}\left({\stackrel{-}{\theta }}_{i+1,j}-{\stackrel{-}{\theta }}_{i-1,j}\right)\right] \\ +\frac{1}{Pr}\left[1+\frac{4}{3N}\right]\right){\stackrel{-}{\mathit{\theta }}}_{\mathit{i}-1,\mathit{j}}+\left(-\frac{{∆Z}^2}{∆X}{X}_i{U}_{i,j}-\frac{2}{Pr}\left[1+\frac{4}{3N}\right]\right){\stackrel{-}{\mathit{\theta }}}_{\mathit{i},\mathit{j}} \\ \hspace{1em}+\left(-\left[\frac{1}{2}∆Z\left(W_{i,j}-\frac{1}{4}{Z}_j{U}_{i,j}\right)-\frac{1}{P_r}\frac{∆Z}{2Z_j}-N_b\frac{1}{4}\left({\stackrel{-}{\phi }}_{i+1,j}-{\stackrel{-}{\phi }}_{i-1,j}\right) \\ \hspace{1em}\hspace{1em} +N_t\frac{1}{4}\left({\stackrel{-}{\theta }}_{i+1,j}-{\stackrel{-}{\theta }}_{i-1,j}\right)\right]+\frac{1}{Pr}\left[1+\frac{4}{3N}\right]\right){\stackrel{-}{\mathit{\theta }}}_{\mathit{i}+,\mathit{j}}=-\frac{{∆Z}^2}{∆X}{X_{i}U}_{i,j}{\stackrel{-}{\theta }}_{i,j-1} \end{gathered}$$

(64)

$$\begin{gathered} \left(\left[\frac{1}{2}∆Z\left(W_{i,j}-\frac{1}{4}{Z}_j{U}_{i,j}\right)+\frac{∆Z}{{2Z}_j}\right]+\frac{1}{L_e}\right){\stackrel{-}{\mathit{\phi }}}_{\mathit{i}-1,\mathit{j}}+\left({-X}_i\frac{{∆Z}^2}{∆X}{U}_{i,j}-\frac{2}{Sc}\right){\stackrel{-}{\mathit{\phi }}}_{\mathit{i},\mathit{j}} \\ +\left(-\left[\frac{1}{2}∆Z\left(W_{i,j}-\frac{1}{4}{Z}_j{U}_{i,j}\right)+\frac{∆Z}{{2Z}_j}\right]+\frac{1}{Sc}\right){\stackrel{-}{\mathit{\phi }}}_{\mathit{i}+{\displaystyle 1},\mathit{j}} \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em} =X_i{U}_{i,j}\frac{{∆Z}^2}{∆X}{\stackrel{-}{\phi }}_{i,j-1}-\frac{1}{Sc}\frac{N_t}{N_b}\left\{\left(1+\frac{∆Z}{{2Z}_j}\right){\stackrel{-}{\theta }}_{i+1,j}+\left(1-\frac{∆Z}{{2Z}_j}\right){\stackrel{-}{\theta }}_{i-1,j}-2{\stackrel{-}{\theta }}_{i,j}\right\} \end{gathered}$$

(65)

Boundary conditions of the modeled problem

$$\begin{gathered} W_{i,j}=U_{i,j}=0,\hspace{1em}{\stackrel{-}{\varphi }}_{i,j}=1,\hspace{1em}{\stackrel{-}{\theta }}_{i,j}=1\hspace{1em}at{Z}_j=0 \\ {U}_{i,j}\to 0,\hspace{1em}{\stackrel{-}{\varphi }}_{i,j}\to 0,{\stackrel{-}{\theta }}_{i,j}\to 0\hspace{1em}as{Z}_j\to \infty \end{gathered}$$

(66)

The approximate solutions determined with the finite difference method are discussed in detail in the forthcoming section. The convergence and mesh sensitivity analysis criterion are presented as follows to achieve accurate numerical solutions for the variables $U,W$, $\mathrm{\theta }$ and $\varphi$, respectively:

$$\max \left|U_{ij}\right|+\max \left|W_{ij}\right|+\max \left|{\mathsf{\theta }}_{ij}\right|+\max \left|{\mathsf{\varphi }}_{ij}\right|\le ϵ$$

where $ϵ={10}^{-5}$. The computation is started at $X=0$ and then marches downstream implicitly. Here, we have taken the step size $\mathrm{\Delta }X=0.05$ and $\mathrm{\Delta }Z=0.02$, respectively.

6. Results and Discussion

The current section is dedicated to the results and discussion in relation to the flow model. The obtained numerical solutions issued from solving the governing equations, around the sphere in region-I and plume region-III, are displayed graphically and in tabular forms in subsequent subsections.

6.1. Analysis of Heat and Fluid Flow Characteristics around the Sphere

The solution of the problem is obtained by solving Equations (1)–(4) with their boundary conditions, as given in Equation (5). A finite difference technique has been applied in order to reach the solution. The study’s primary goal is to examine the effects of various flow model inputs, on the velocity profile $U$, temperature distribution $\theta$, and nanoparticle concentration $\varphi$. In addition, the quantitative results are presented in terms of skin friction coefficient $C_{f}G{r}^{1/4}$, Nusselt $NuG{r}^{-1/4}$, and the Sherwood number $ShG{r}^{-1/4}$. The pertinent parameters of interest are the radiation parameter $N$, the thermophoresis parameter $N_t$, the Brownian motion $N_b$, the Schmidt number Sc, and the Prandtl number Pr. Figure 2, Figure 3 and Figure 4 show the effect of varying the thermal radiation parameter $N$ on velocity profile $U$, temperature distribution $\theta$, and nanoparticle concentration $\varphi$, respectively. It is noticed that as $N$ is raised, velocity, temperature, and nanoparticle concentrations increase for all the considered positions around the surface of a sphere. The effects of N are included in order to increase the fluid flow domain temperature, as shown in the aforementioned graphs. Figure 2, Figure 3 and Figure 4 illustrate how an increase in N raises the fluid’s thermal conductivity and mean absorption efficiency. Both contribute to raising the temperature of the fluid flow domain. As the fluid’s temperature increases, the fluid thins out and the resistance between the layers increases. As a result, the velocity field is increased in accordance. It is worthy to note that the velocity profile $U$ is largest at point $X$ = 2.0 and that the velocity and concentration profiles are maximum at position $X$ = 3.0. The overall established trend is that the increase in N leads to an increase in $\theta$ and $\varphi$. Figure 5, Figure 6 and Figure 7 demonstrate the solutions for $U$, $\theta$ and $\varphi$ profiles, respectively, for different values of thermophoresis parameter $N_t$ in the sphere region. It is observed that all the variables become stronger when there is an increase in the magnitude of $N_t$. In these figures, we can conclude that $U$ is maximum at position X = 2.0, while $\theta$ and $\phi$ are maximum at position X = 3.0 in the sphere region. The observed trends are physically true because the increment in $N_t$ is caused by the increase in temperature difference between the surface temperature and ambient temperature; the viscous force between the layers of the fluid weakens, causing velocity to become boosted and the temperature of the fluid to be enhanced. Figure 8, Figure 9 and Figure 10 show how the Brownian motion parameter $N_b$ affects velocity, temperature distribution, and nanoparticle concentration, respectively. All the variables increase against the increasing values of $N_b$. This result is justified by the fact that an increment in $N_t$ is caused by an increase in nanoparticle concentration difference between the surface and ambient concentration, and viscous force is rescued, causing the increment in concentration. It can also be noticed that the contribution of $Nb$ is maximum at point X = 2.0 for $U$ and at X = 3.0 for $\theta$ and $\phi$ for Nb = 1.5. The effects of Sc on dynamical and thermal behaviors are presented in Figure 11, considering that the other parameters are constant. From the analysis of these figures, it can be concluded that $\varphi$ increases due to the increase in the Brownian motion, and it reaches its maximum at point X = 0.1.

Table 1. Outcomes of the Schmidt number $Sc$ on (a) the skin friction coefficient $C_{f}G{r}^{1/4}$, (b) Nusselt $NuG{r}^{-1/4}$, and (c) the Sherwood number $ShG{r}^{-1/4}$ when $Nt=0.2,Nb=10.0,N=1.1,$ and $Pr=7.0.$

X	${\mathit{C}}_{\mathit{f}}\mathit{G}{\mathit{r}}^{1/4}$		$\mathit{N}\mathit{u}\mathit{G}{\mathit{r}}^{-1/4}$		$\mathit{S}\mathit{h}\mathit{G}{\mathit{r}}^{-1/4}$
	$\mathit{S}\mathit{c}=0.1$	$\mathit{S}\mathit{c}=2.0$	$\mathit{S}\mathit{c}=0.1$	$\mathit{S}\mathit{c}=2.0$	$\mathit{S}\mathit{c}=0.1$	$\mathit{S}\mathit{c}=2.0$
0.1	0.00981	0.01145	0.06308	0.06349	0.07619	0.07801
1.0	0.06199	0.04599	0.06578	0.08904	0.07682	0.10350
2.0	0.10721	0.04818	0.19978	0.09024	0.20579	0.10476
3.1	0.00534	0.00578	0.05831	0.05809	0.07276	0.07327

,

Table 2. Outcomes of the Prandtl number $Pr$ on (a) the skin friction coefficient $C_{f}G{r}^{1/4}$, (b) Nusselt $NuG{r}^{-1/4}$, and (c) the Sherwood number $ShG{r}^{-1/4}$ when $Nt=0.2,Nb=10.0,N=1.1,$ and $Sc=10.0.$

X	${\mathit{C}}_{\mathit{f}}\mathit{G}{\mathit{r}}^{1/4}$		$\mathit{N}\mathit{u}\mathit{G}{\mathit{r}}^{-1/4}$		$\mathit{S}\mathit{h}\mathit{G}{\mathit{r}}^{-1/4}$
	$\mathit{P}\mathit{r}=7.0$	$\mathit{P}\mathit{r}=10.0$	$\mathit{P}\mathit{r}=7.0$	$\mathit{P}\mathit{r}=10.0$	$\mathit{P}\mathit{r}=7.0$	$\mathit{P}\mathit{r}=10.0$
0.1	0.01271	0.01483	0.06050	0.05501	0.07456	0.07243
1.0	0.08087	0.08239	0.10321	0.09503	0.12380	0.12105
2.0	0.08619	0.08762	0.10548	0.09702	0.12654	0.12361
3.1	0.00611	0.00742	0.05663	0.04794	0.07183	0.06604

,

Table 3. Outcomes of the radiation parameter $N$ on (a) the skin friction coefficient $C_{f}G{r}^{1/4}$, (b) Nusselt $NuG{r}^{-1/4}$, and (c) the Sherwood number $ShG{r}^{-1/4}$ when $Nt=0.2,Nb=0.3,Sc=10.0,$ and $Pr=7.0.$

X	${\mathit{C}}_{\mathit{f}}\mathit{G}{\mathit{r}}^{1/4}$		$\mathit{N}\mathit{u}\mathit{G}{\mathit{r}}^{-1/4}$		$\mathit{S}\mathit{h}\mathit{G}{\mathit{r}}^{-1/4}$
	$\mathit{N}=0.7$	$\mathit{N}=1.1$	$\mathit{N}=0.7$	$\mathit{N}=1.1$	$\mathit{N}=0.7$	$\mathit{N}=1.1$
0.1	0.01812	0.01421	0.01271	0.06050	0.07342	0.07456
1.0	0.15241	0.08169	0.08087	0.10321	0.12246	0.12380
2.0	0.16466	0.08695	0.08619	0.10548	0.12508	0.12654
3.1	0.02562	0.00611	0.05072	0.05663	0.00679	0.07183

,

Table 4. Outcomes of the thermophoresis parameter $N_t$ on (a) the skin friction coefficient $C_{f}G{r}^{1/4}$, (b) Nusselt $NuG{r}^{-1/4}$, and (c) the Sherwood number $ShG{r}^{-1/4}$ when $Sc=10.0,Nb=0.3,N=1.0$, and Pr = 7.0.

X	${\mathit{C}}_{\mathit{f}}\mathit{G}{\mathit{r}}^{1/4}$		$\mathit{N}\mathit{u}\mathit{G}{\mathit{r}}^{-1/4}$		$\mathit{S}\mathit{h}\mathit{G}{\mathit{r}}^{-1/4}$
	$\mathit{N}\mathit{t}=0.2$	$\mathit{N}\mathit{t}=0.9$	$\mathit{N}\mathit{t}=0.2$	$\mathit{N}\mathit{t}=0.9$	$\mathit{N}\mathit{t}=0.2$	$\mathit{N}\mathit{t}=0.9$
0.1	0.01271	0.01093	0.06050	0.00835	0.07456	0.02171
1.0	0.08087	0.25808	0.10321	2.87762	0.12380	0.18380
2.0	0.08619	0.26958	0.10548	0.09417	0.12654	0.19084
3.1	0.00611	0.01412	0.05663	0.06653	0.07183	0.05583

and

Table 5. Outcomes of the Brownian motion parameter $N_b$ on (a) the skin friction coefficient $C_{f}G{r}^{1/4}$, (b) Nusselt $NuG{r}^{-1/4}$, and (c) the Sherwood number $ShG{r}^{-1/4}$ when $Nt=0.2,Sc=10.0,N=1.1,$ and $Pr=7.0$.

X	${\mathit{C}}_{\mathit{f}}\mathit{G}{\mathit{r}}^{1/4}$		$\mathit{N}\mathit{u}\mathit{G}{\mathit{r}}^{-1/4}$		$\mathit{S}\mathit{h}\mathit{G}{\mathit{r}}^{-1/4}$
	$\mathit{N}\mathit{b}=0.3$	$\mathit{N}\mathit{b}=0.9$	$\mathit{N}\mathit{b}=0.3$	$\mathit{N}\mathit{b}=0.9$	$\mathit{N}\mathit{b}=0.3$	$\mathit{N}\mathit{b}=0.9$
0.1	0.01271	0.03029	0.06050	0.03671	0.07456	0.07692
1.0	0.08087	0.12463	0.10321	0.06132	0.12380	0.11358
2.0	0.08619	0.13128	0.10548	0.06242	0.12654	0.11540
3.1	0.00611	0.01716	0.05663	0.02969	0.07183	0.07301

, illustrate the variations in the physical quantities: skin friction $C_{f}G{r}^{1/4}$_, the Nusselt number $NuG{r}^{-1/4}$, and the Sherwood number $ShG{r}^{-1/4}$ for diverse values of the governing parameters. These results are evaluated at various locations on the surface of the sphere. Table 1 presents the results of $C_{f}G{r}^{1/4}$, $NuG{r}^{-1/4}$, and $ShG{r}^{-1/4}$ for various values of the Schmidt number $Sc$ at different points of the surface of the sphere. It is noticed that $Sc$, $C_{f}G{r}^{1/4}$, $NuG{r}^{-1/4}$, and $ShG{r}^{-1/4}$ are higher at $X=\mathrm{0.1,3.1}$ and lower at $X=1.0,2.0$. On the other hand, $NuG{r}^{-1/4}$ and $ShG{r}^{-1/4}$ risefor $X=\mathrm{0.1,2.0,3.1}$ but reduceat $X=1.0$. Similarly, Table 2 depicts the effect of the Prandtl number on the evaluated physical properties. It has to be mentioned that the increase in $Pr$ causes the increase in $C_{f}G{r}^{1/4}$, but reductions occur in $NuG{r}^{-1/4}$ and $ShG{r}^{-1/4}$. In Table 3, the effects of $N$ on $C_{f}G{r}^{1/4}$, $NuG{r}^{-1/4}$, and $ShG{r}^{-1/4}$ are illustrated. Augmenting $N$ causes a decline in $C_{f}G{r}^{1/4}$ and an increase in $NuG{r}^{-1/4}$ and $ShG{r}^{-1/4}$. Table 4 portrays the trend of the same properties for different values of $N_t$. It can be seen that all the properties go down by increasing the values of $N_t$. Table 5 presents the behavior of the skin of friction coefficient, Nusselt number, and Sherwood number for various values of $N_b$. It is observed that each variable of interest increases against increasing values of $N_b$.

6.2. Physical Behavior of Material Properties in the Plume Region-III

The objective of the present work is to investigate the nanofluid boundary layer flow around the sphere and the eruption of the fluid into the plume region-III, as seen in Figure 1. For this purpose, some assumptions are considered in the plume region, and the issued mathematical model is given in Equations (34)–(37) along with the specified boundary conditions (Equation (38)). After making necessary arrangements regarding the scaling of the variables, we transformed the coupled system of the partial differential equation into a convenient form for the integration. A numerical analysis of the thermofluid dynamics and nanoparticle motion characteristics is conducted by adopting the finite difference method, along with the Gaussian elimination algorithm to calculate the values of the unknown variables $U$,$\theta$, and $\mathrm{\varphi }$ for different values of the radiation parameter $N$. In Figure 12, Figure 13 and Figure 14, the effects of radiation on the fluid dynamics, thermal field, and nanoparticle motion are displayed. Fluid velocity $U$, thermal field $\theta$, and nanoparticle motion $\varphi$ decreased drastically as the radiation parameter $N$ increased. This is because with the increase in $N$, the fluid becomes grossly optically dense gray. On the other hand, the nanoparticle motion restored after retardation.

Table 6. Outcomes of radiation parameter $N$ on (a) the skin friction coefficient $C_{f}G{r}^{1/4}$, (b) Nusselt $NuG{r}^{-1/4}$, and (c) the Sherwood number $ShG{r}^{-1/4}$ when $Nt=0.2$, $Nb=0.3,Sc=1.0$, and $Pr=7.0$.

X	${\mathit{C}}_{\mathit{f}}\mathit{G}{\mathit{r}}^{1/4}$		$\mathit{N}\mathit{u}\mathit{G}{\mathit{r}}^{-1/4}$		$\mathit{S}\mathit{h}\mathit{G}{\mathit{r}}^{-1/4}$
	$\mathit{N}=0.5$	$\mathit{N}=0.7$	$\mathit{N}=0.5$	$\mathit{N}=0.7$	$\mathit{N}=0.5$	$\mathit{N}=0.7$
0.0	50.75594	46.76019	0.95139	0.23755	192.65468	179.60024
0.1	32.08902	34.08479	7.68669	7.03547	82.65567	75.70414
1.0	15.67088	13.73597	9.27360	11.81349	0.11835	0.12578
2.0	10.61131	8.79781	9.27483	11.81388	0.11340	0.11728
3.0	7.79274	5.95524	18.88779	17.98556	3.34746	4.44476
4.0	5.76776	3.84707	18.63854	18.20510	2.30285	4.74118
5.0	4.13138	2.18996	18.63854	17.75627	2.29610	6.89948
6.0	2.79540	1.27903	17.91728	16.67107	3.51993	9.09658
7.0	1.68713	0.95814	16.81319	15.59305	6.43594	12.18080
8.0	1.14654	0.78669	16.92500	14.97582	9.41521	13.28487
9.0	0.91312	0.68151	15.51050	14.72355	11.44516	13.84803
10.0	0.77717	0.60034	15.29638	14.52141	12.45936	14.1815

describes the behavior of different scaling parameters on the physical quantities of the skin friction coefficient $C_{f}G{r}^{1/4}$, Nusselt $NuG{r}^{-1/4}$, and the Sherwood number $ShG{r}^{-1/4}$ in plume region-III at the leading edge. We can see that all the properties decline when the value of $N$ increases.

7. Conclusions

The current study analyzes the flow of nanofluid in a laminar boundary layer around a heated sphere inside a plume while thermal radiation is present. It can be concluded that the temperature distribution and mass concentration were improved when $N,N_t,Sc,$ and $N_b$ were increased. It is also noticed that as $Sc,Pr$ and $N_b$ were increased, the skin friction became stronger, but an opposite trend was found for the increasing values of $N$ and $N_b$. Tabular results showed that increasing the values of $Sc,Pr$ and $N$ resulted in an enhancement of the heat transfer rate. An opposite behavior was encountered when $Nt$ and $Nb$ were increased. The numerical results that corresponded with the nanoparticle concentration showed that as $Sc$ and $Nb$ were increased, an increasing trend in nanoparticle concentration was noticed. It is also interesting to mention that the increase in $N$ and $N_t$ led to a reduction in the above-discussed properties. In the plume region-III, it was noticed that increasing the values of $N$ caused a reducing trend in velocity and temperature fields, thus leading to a rising trend for the concentration. A decreasing trend was also found to occur for $C_{f}G{r}^{1/4}$, $NuG{r}^{-1/4}$, and the Sherwood number $ShG{r}^{-1/4}$.