# Fully Developed Flow of a Nanofluid through a Circular Micropipe in the Presence of Electroosmotic Effects

Authors: Muhammad D. K. Niazi; Hang Xu
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

Microscale heat sinks based on channels or pipes are designed to restrict the temperatures of microelectromechanical systems, which have a wide range of applications in the modern engineering and mechanics. In this context, this work aims to study heat convection and entropy generation of a fully developed nanofluid flow in a circular micropipe in the presence of an electrical double layer. Buongiorno’s model is employed to exhibit the nanofluid behavior. The governing equations are reduced to a system of nonlinear ordinary differential equations through appropriate similarity transformations. Particularly, we rectify the pressure term as an unknown constant, which makes our flow model compatible with those well-known fluid flow models in macrosize. Highly accurate solutions are obtained and verified. Analysis for physical properties of electric field, velocity field, temperature, and nanoparticle distributions is discussed followed by an investigation of the entropy evolution in the flow. The results show that flow behavior and total entropy of the system depend on the electroosmosis, thermophoresis, and fluid viscosity. However, the influence of the electrical double layer on the flow and system entropy is negligible when the electroosmotic parameter exceeds a maximum value.


---

## Body

## 1. Introduction

Microelectromechanical systems (MEMS) are of great significance due to their vast applications in modern electronic equipment. In some of these devices, order of heat dissipation is expected to be as high as100W/cm2 [1], which could cause breakdowns and failures. In this context, cooling pumps with microsized channels or pipes are deployed. Tuckerman and Pease [2] were the first to suggest the idea of microchannel flow as heat sinks. Their concept was then adopted by several researchers [3, 4] to study the thermal properties and behaviors in such channels/pipes. Particularly, Phillips [5] concluded that the thermal conductivity and the viscosity are two essential properties that affect heat transfer behavior of microchannel flows.Nanoscale particles, with a diameter of 1 nm–100 nm, are classified as nanoparticles, and a combination of these particles with a base liquid is known as nanofluid. Choi and Eastman [6] were the first to suggest the use of nanoparticles with the base fluid to improve its thermal conductivity. Their approach has been used to enhance the thermal conductivity of fluids in various studies [7–10]. Specifically, a wide range of particles have been used to observe their specific behavior in a base fluid [11, 12]. From a theoretical point of view, several researchers have attempted to mathematically model the nanoparticle behavior. Particularly, Jang and Choi [13] concluded the existence of a link between nanoparticle diffusion and thermal properties of the nanofluid and proposed a thermal conductivity model based on temperature. Various other models were designed to study nanoparticle behavior in nanofluid, where the two-phase models [14–16] and the Tiwari–Das model [17] have been applied widely in the literature.Buongiorno [18] proposed a better approach in 2006 which incorporated nanoparticle and base fluid properties to obtain a more accurate model. In his study, the fluid layers in the neighborhood of a nanoparticle were considered as a continuum and seven mechanisms of molecular disposition were used. However, thermophoretic and Brownian diffusion are two most important phenomena considered in his model [19]. In contrast to Buongiorno’s model, the Tiwari–Das model uses nanoparticle volumetric fraction, whereas the two-phase models use physical phenomena such as gravity, interphase friction, sedimentation, Brownian diffusion, and dispersion. However, these approaches do not consider sufficient effects such as thermophoresis [20]. In this context, Tsai [21] considers mechanisms such as electrophoresis and convection in addition to Brownian and thermophoretic diffusion as the major primary nanoparticle deposition mechanisms. Therefore, Buongiorno’s model has been used in various studies to model nanoparticle behavior. In particular, Kuznetsov and Nield [22] studied the nanoparticles’ impact on a natural convection flow on a vertical plate, while Hwang et al. [23] studied their impact in a laminar forced-convection flow. Nield and Kuznetsov [24, 25], Malvandi et al. [26], Sheikhzadeh et al. [27], and Saleem et al. [28] performed extensive studies on nanofluid convective transport. All these studies highlight the adequacy of Buongiorno’s model.Mala et al. [29] experimentally confirmed the significance of the electrokinetic effects and other factors such as channel material and size on the characteristics of microchannel flow. Electroosmosis is an electrokinetic effect defined as the motion of the ionized liquid in relation to stationary charged surfaces, and one of the most important mechanisms due to this motion is the electrical double layer. Several theoretical and experimental studies investigating the electroosmotic effects on flow behavior can be found in the literature [30, 31]. In particular, Hunter [32] observed that electrical double layer (EDL) near the microchannel walls are generated by mutual attraction between free ions in the fluid and the static ions in the channel boundaries. Since most of the MEMS operate in the presence of electrical double layer, it is essential to investigate their influence on flow and heat transfer behavior. Mala et al. [33] modeled this behavior with the assumption that the electrostatic potential vanishes in the middle part of the channel. However, this assumption is only valid when the height between channel walls is very large. Ren et al. [34] made an improvement to this model with the symmetric flow pattern being assumed. From then on, many researchers [35, 36] considered various aspects of microchannel flow behavior.In theoretical and experimental studies on microchannels, the pressure term was commonly prescribed as a constant. This is true in real circumstance. However, this term is closely associated with other physical quantities such as velocity, temperature, and fluid properties. Once it is prescribed, several other parameters are automatically determined. For example, in the well-known Poiseuille flow problem, the constant pressure is closely associated with the Reynolds number; we can determine the flow pattern either by the constant pressure or by the Reynolds number when the fluid properties remain unchanged. However, most researchers were not aware of this issue, which caused their models to be incompatible. Recently, Xu et al. [37] overcame this incompatibility in their model by assuming that the Reynolds number is known while the pressure constant is unknown. Their proposed model was further testified in a series of alternative microchannel configurations by Xu and his collaborators [38, 39].There has been growing interest in optimization of system entropy in order to achieve better performance. According to the second law of thermodynamics, every real process is irreversible. However, this amount of irreversibility can be evaluated through entropy generation rate of the system. It is known that the quality of energy is lower for a higher entropy system [40]. Bejan [41] suggested that heat transfer and viscous friction are the major factors of entropy generation in a forced-convection viscous flow. Many researchers [42, 43] made efforts to understand the mechanism of entropy generation. Apart from the factors of heat transfer and viscous friction, other parameters such as Joule heating, electromagnetic diffusion, magnetic field, and nanoparticle distribution could be significant to affect energy generation in a microchannel system. The corresponding analysis was performed by Xie and Jian [44] who considered six factors that include heat transfer, viscous friction, Joule heating, electromagnetic diffusion, magnetic field, and nanoparticle distribution to analyze the fluid irreversibility rate in a microchannel with both magnetic and electrical fields being taken into account simultaneously. Xie and Jian considered a two-layer magnetohydrodynamic electroosmotic microchannel flow and analyzed the entropy of the system. A similar research was conducted by Escandón et al. [45], who studied entropy generation of non-Newtonian fluids affected by electroosmotic effects in a microchannel. Earlier in 2010, Zhao and Liu [46] studied evolution of entropy generation in open-end and closed-end microchannel flows under influence of electroosmosis. Some recent studies on entropy analysis are from Hayat et al. [47], Rashid et al. [48], and Zeeshan et al. [49].Motivated from the aforementioned research, this work is an investigation of entropy generation in a fully developed nanofluid flow in a microradius circular pipe in the presence of electrical effects from an EDL. To the best of the authors’ knowledge, entropy analyses of electroosmotic nanofluid flow in a microradius pipe are rare in the literature. However, this research is necessary due to wide range of applications of MEMS and the importance of microradius structures in the design of comparable heat sinks. Buongiorno’s nanofluid model is used to describe the nanofluid flow behavior. The pressure term is assumed as an unknown constant as proposed by Xu et al. [37]. The governing equations embodying the conservation of electrical potential, total mass, momentum, energy, and nanoparticle volume fraction are reduced to a system of ordinary differential equations through appropriate similarity transformations. Highly accurate solutions are obtained by the improved homotopy analysis method (HAM) [50]. The physical behaviors of the flow are presented and discussed. The dimensionless parameters’ effects on the system entropy are analyzed in order to optimize system performance.


## 2. Mathematical Modeling

Consider a laminar, fully developed flow of a nanofluid driven in a microradius pipe by a constant pressure gradient in the presence of an electrical double layer. The pipe is placed along thez-axis with its radius being R and its length being L. It is assumed that the length is large enough so that the effects from the openings at the end of the pipe can be eluded. The physical sketch is presented in Figure 1 with the cylindrical coordinate system r,ϕ,z being used.Figure 1
Physical sketch of the problem.Because of the small radius of the pipe, parallel flow hypothesis can be assumed. The flow is assumed to be symmetrical in angular direction, i.e.,∂/∂ϕ=0. Since the zeta potential is constant throughout the pipe surface, the electrical potential depends only on the radial component, i.e., ψ=ψr (refer to [51]). It shall be noted that the velocities in r-, ϕ-directions, based on the parallel flow hypothesis [52] (i.e., the length L is far greater than the radius R of the pipe), are negligible, i.e., u≈0 and v≈0. Hence, only the z-direction velocity component w is considered. The flow is thermally fully developed with a constant wall heat flux which implies that the temperature gradient in z-direction is constant, i.e., ∂T/∂z=C1 (refer to [53]). Furthermore, the nanoparticle mass flux on the wall is assumed to be constant, i.e., ∂C/∂z=C2. The induced electric field vector is E=0,0,Ez, while the flow is considered in the absence of external electric and magnetic fields. Since the flow is laminar, inertia effect, diffusiophoresis, Magnus effect, fluid drainage, and gravity effect can be neglected when considering Buongiorno’s model [18].With those assumptions, the continuity equation is automatically satisfied. Hence, the governing equations embodying the conservation of the electrostatic potential, total mass, momentum, energy, and nanoparticle volume fraction are written as(1)1r∂∂rr∂ψ∂r=−ρeε,(2)−∂p∂z+μr∂∂rr∂w∂r+Ezρe=0,(3)w∂T∂z=α1r∂∂rr∂T∂r+μρcf∂w∂r2+σEz2ρcf+τDB∂T∂r∂C∂r+∂T∂z∂C∂z+τDTT0∂T∂r2+∂T∂z2,(4)DB1r∂∂rr∂C∂r+DTT01r∂∂rr∂T∂r−w∂C∂z=0,where ψ is the electric potential, w is the velocity in z-direction, T is the temperature, C is the nanoparticle volume fraction, Ez is the streaming potential in z-direction, α=kf/ρcf is the thermal diffusivity, τ=ρcs/ρcf is the heat capacity ratio, ρ is the fluid density, c is the specific heat, kf is the thermal conductivity with the subscripts f and p denoting the quantity for the fluid or solid particles, respectively, ρe is the charge density, ε is the fluid medium permittivity, p is the pressure, μ is the viscosity, σ is the electrical conductivity, T0 is the reference temperature, DB is the coefficient of Brownian diffusion, and DT is the coefficient of thermophoretic diffusion.The flow is assumed to satisfy the no-slip condition, and the zeta potential, the temperature, and the nanoparticle volume fraction on the boundary are assumed to be constant, given asζ∗, Tw, and Cw, respectively. So, the appropriate boundary conditions for equations (1)–(4) are(5)r=R:ψ=ζ∗,w=0,T=Tw,C=Cw,r=0:dψdr=0,∂w∂r=0,dTdr=0,dCdr=0,along with the average velocity condition due to the conservation of flow flux:(6)Wm=2R∫0Rurdr.With a uniform dielectric constant and the assumption that the thermal energy is significantly larger than the electrical potential, the Poisson–Boltzmann equation (1) is reduced to [54](7)1r∂∂rr∂ψ∂r=k2ψ,where k2=2n0z^2e2/εkbT^ is the Debye–Hückel parameter, whose reciprocal determines the width of EDL. Here, n0 is the bulk concentration of ions, z^ is the ion valence, e is the fundamental charge, kb is the Boltzmann constant, and T^ is the absolute temperature.


## 3. Nondimensional Reduction

The following similarity variables are defined to reduce the governing equations to a nondimensional form:(8)η=rR,κ=Rk,θ=T−T0Tw−T0,Ψ=z^eψkbT^,W=wWm,s=C−C0Cw−C0.It shall be noted that the domain [0, R] is transformed to [0, 1] after application of transformations (8), and the respective governing equations become(9)Ψ″+Ψ′η−κ2Ψ=0,(10)W″+W′η+Γ−2G1E¯sΨ=0,(11)θ″+θ′η+Ntθ′2+Nbs′θ′+BrW′2+G2E¯s2+A1W+G3NbA2+G4NtA1=0,(12)s″+s′η+NtNbθ″+θ′η+A2W=0,which are subject to the following boundary conditions:(13)r=1:Ψ=ζ,W=0,θ=1,s=1,(14)r=0:Ψ′=0,W′=0,θ′=0,s′=0,and the average velocity condition becomes(15)2R∫01Wηdη=1.Here,ζ=z^eζ∗/kbT^ represents the dimensionless zeta potential and κ, Γ, E¯s, Nt, Nb, Br, A1, and A2 denote the nondimensional parameters for the electroosmosis, i.e., the pressure gradient, the streaming potential, the thermophoretic diffusion, the Brownian diffusion, the Brinkman number, the constant temperature, and the nanoparticle concentration gradients in z-direction, respectively, which are defined as(16)κ=Rz^e2n0εkbT^,Γ=−R2p′zμWm,Nt=τDTTw−T0T0α,Nb=τDBCw−C0α,Br=μWm2kfTw−T0,A1=−R2WmC1αTw−T0,A2=−R2WmC2DBCw−C0,G1=ζR2n0z^eμLWm,G2=σR2ζ2L2kfTw−T0,G3=DBαR2Wm2,G4=α2R2Wm2.The nondimensional streaming potential parameterE¯s in equation (11) is unknown, which is derived by using the equilibrium condition of the streaming current Is produced from the movement of charged particles in the diffuse layer and the conduction current Ic produced by the streaming potential generated from the streaming current (refer to [51]), which is obtained as(17)E¯s=4G5∫01ηWΨdη,where G5 = LWmn0z^e/(ζλ0).


### 3.1. Entropy Analysis

Entropy of the system can be computed after the solutions for the system of equations (9)–(14) are obtained. According to Bejan’s conclusion [42], the entropy generation rate is given by(18)SGG=SGH+SGJ+SGF+SGC,where SGH, SGJ, SGF, and SGC are the entropy terms for the reversibility due to heat transfer, Joule heating, viscous friction, and the nanoparticle distribution, respectively, which are written as(19)SGH=kfT2∂T∂r2+∂T∂z2,SGJ=σEz2T,SGF=μTΦ,SGC=RDC∂C∂r2+∂C∂z2+RDT∂C∂r∂T∂r+∂C∂z∂T∂z.The similarity variables defined in (8) are used to obtain dimensionless form of equation (18). The characteristic entropy rate SG=R2SGG/kf is given as(20)SG=1θ+θ01θ+θ0θ′2+G4A12+BrW′2+G2E¯s2+MDs+Λ0s′2+G32G4A22+MDθ+θ0s′θ′+G3A1A2.In the above expression,MD, θ0, and Λ0 are given as(21)MD=RDCw−C0kf,θ0=T0Tw−T0,Λ0=C0Cw−C0.Finally, the total entropy in the pipe can be computed by taking an integral over the pipe radius:(22)Stotal=∫01SGηdη.


## 3.1. Entropy Analysis

Entropy of the system can be computed after the solutions for the system of equations (9)–(14) are obtained. According to Bejan’s conclusion [42], the entropy generation rate is given by(18)SGG=SGH+SGJ+SGF+SGC,where SGH, SGJ, SGF, and SGC are the entropy terms for the reversibility due to heat transfer, Joule heating, viscous friction, and the nanoparticle distribution, respectively, which are written as(19)SGH=kfT2∂T∂r2+∂T∂z2,SGJ=σEz2T,SGF=μTΦ,SGC=RDC∂C∂r2+∂C∂z2+RDT∂C∂r∂T∂r+∂C∂z∂T∂z.The similarity variables defined in (8) are used to obtain dimensionless form of equation (18). The characteristic entropy rate SG=R2SGG/kf is given as(20)SG=1θ+θ01θ+θ0θ′2+G4A12+BrW′2+G2E¯s2+MDs+Λ0s′2+G32G4A22+MDθ+θ0s′θ′+G3A1A2.In the above expression,MD, θ0, and Λ0 are given as(21)MD=RDCw−C0kf,θ0=T0Tw−T0,Λ0=C0Cw−C0.Finally, the total entropy in the pipe can be computed by taking an integral over the pipe radius:(22)Stotal=∫01SGηdη.


## 4. Solution Procedure

Solutions of equations (9)–(14) are obtained using the homotopy analysis method (HAM) with the help of Mathematica computational software. The first step of solving a system of differential equations by the HAM is the formulation of the deformation equation which is given as(23)1−qℒfFη;q−f0η=qℏFNfFη;q,where Fη;q represents a mapping function of fη, ℒF denotes the chosen linear operator, NF denotes the nonlinear operator, the embedding parameter q∈0,1, and ℏF denotes the convergence control parameter. From the deformation equation, the mth-order HAM deformation equation can be obtained by taking m-times derivative with respect to q and dividing the resulting equation by m!, then finally letting q=0. So, from equation (23), it becomes(24)ℒFfmη−χmf0η=ℏFRmη,where(25)fmη=1m!∂mFη;q∂qmq⟶0,Rmη=1m−1!∂m−1NFη;q∂qm−1q⟶0,χm=0,m≤1,1,m>1.In order to solve the problem in the HAM framework, the unknown solutions forΨ, W, θ, and s are assumed to be of the following form:(26)Ψ=Ψ0+∑j=1∞Ψj,W=W0+∑j=1∞Wj,θ=θ0+∑j=1∞θj,s=s0+∑j=1∞sj,and the unknown pressure gradient constant Γ is assumed to be(27)Γ=Γ0+∑j=1∞Γj.The linear operator is chosen to be(28)ℒΨ=ℒW=ℒθ=ℒs=ℒ=∂2∂η2.With the assumed solution series expressions (26) and the linear operator (28), the mth-order HAM deformation equations can be written as(29)Ψm″−χmΨm−1″=ℏΨRΨ,m,Wm″−χmWm−1″=ℏWRW,m,θm″−χmθm−1″=ℏθRθ,m,sm″−χmsm−1″=ℏsRs,m,subject to the boundary conditions(30)η=0:Ψm′=Wm′=θm′=sm′=0,η=1:ϕm=Wm=θm=sm=0,with the following additional condition:(31)∫01ηWmηdη=0.The initial guesses are chosen as(32)Ψ0η=ζ,θ0η=s0η=1,W0η=2R1−η2.In equation (29), ℏΨ, ℏW, ℏθ, and ℏs denote the convergence control parameters, whereas RΨ,m, RW,m, Rθ,m, and Rs,m are(33)RΨ,m=Ψm−1″−κ12Ψm−1+Ψm−1′η,RW,m=Wm−1″+Wm−1′η+Γm−1−8G1G5∑i=0m−1ΩiΨm−i−1,Rθ,m=θm−1″+θm−1′η+∑i=0m−1Nbθi′sm−1−i′+Ntθi′θm−1−i′+BrWi′Wm−1−i′+16G2G5∑i=0m−1Ωi′Ωm−1−i′+A1Wm−1+1−χmG3NbA2+G4NtA1,Rs,m=sm−1″+sm−1′η+NtNbθm−1″+θm−1′η+A2Wm−1,where Ωi=∫01η∑j=0iWjΨi−jdη. The solutions of equation (29) are expressed as follows:(34)Ψmη=Ψ∗η+χmΨm−1η+C1,m+C2,mη,(35)Wmη=W∗η+χmWm−1η+C3,m+C4,mη,(36)θmη=θ∗η+χmθm−1η+C5,m+C6,mη,(37)smη=s∗η+χmsm−1η+C7,m+C8,mη,where Ψm∗η, Wm∗η, θm∗η, and sm∗η are particular solutions, defined by(38)Ψ∗=ℒ−1ℏΨRΨ,m,W∗=ℒ−1ℏWRW,m,θ∗=ℒ−1ℏθRθ,m,s∗=ℒ−1ℏsRs,m,and ℒ−1 is the inverse linear operator which holds the following properties:(39)ℒ−1c0ηm=c0ηmm+2m+1,ℒ−1a0η+b0η=ℒ−1a0η+ℒ−1b0η.The constantsCi,m, i=1,2,…,8, defined in equations (34)–(37) can be obtained using the boundary conditions (30). It shall be noted that the streaming potential parameter E¯s and the pressure gradient constant Γ are unknowns in equations (10) and (11) which are computed during the HAM computations. In order to compute Γ, the additional condition (15) is used on each iteration, and as a result, the series approximation for Γ is obtained. E¯s is obtained in a similar manner, but since the series approximation is not needed, the solution from the highest order computation is used in equation (17).


## 5. Results

This section is divided into two parts. In the first part, the HAM solution is verified using the respective errors from the equations and comparing the electric potential obtained through the analytical and HAM approach. The latter part is a discussion on the physical aspects of the obtained solutions.


### 5.1. Solution Validation

The convergence of the HAM technique has been tested and is presented in Figure2. The 20th-order HAM solutions are used to determine the convergence domains for the parameters. It shall be noted that, for every convergence control parameter, there exists a convergence domain. Therefore, it can be stated that the applied technique is convergent. For the validation of the obtained solution, the approach used by [54] is used. Accuracy of the solution is checked using the maximum total average squared error which is given as(40)Em=maxEΨm,EWm,Eθm,Esm,where m represents the order of the HAM computation and EΨm,EWm,Eθm,andEsm denote the functions for error evaluation, which are defined, based on equations (9)–(12), given as(41)EΨm=∫01Ψ″+Ψ′η−κ2Ψ2dη,EWm=∫01W″+W′η+Γ−2G1E¯sΨ2dη,Eθm=∫01θ″+θ′η+Ntθ′2+Nbs′θ′+BrW′2+G2E¯s2+A1W+G3NbA2+G4NtA1,Esm=∫01s″+s′η+NtNbθ″+θ′η+A2W2dη.Figure 2
HAM convergence control parameter curves. Here,R=1, ζ=1, κ=1, Nb=Nt=1, Br=1, A1=A2=1, and G1=G2=G3=G4=G5=1.Equation (40) is used to compute the mth-order HAM computation error. Table 1 illustrates the errors computed using the following parameter values: R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1. Note in the table that the error decreases as the HAM computational order is increased; this shows that the results are convergent. It should be noted that the homotopy-Padé technique [50] was used to increase the convergence rate and obtain more accuracy, indicating the insignificance of the convergence control parameters which can hence be given any nonzero value. The 50th-order HAM solutions have a good enough accuracy to be used for the foregoing analysis. The HAM approach for solution computation can be further verified using the analytical solution obtained for equation (9). After application of the boundary conditions (15) and (16) and neglecting the imaginary parts, the solution is obtained as(42)ψη=ζI0κrI0κ,where I0 is the zeroth-order modified Bessel function of the first kind.Table 1
The maximum errorEm.
Order
κ
=
1
κ
=
3
κ
=
5
κ
=
10
10
41.766
8.94
×
10
6
2.08
×
10
6
794.505
20
0.52
0.025
8.063
3.651
30
5.138
×
10
−
7
4.803
×
10
−
7
8.2
×
10
−
4
2.557
40
5.877
×
10
−
12
5.552
×
10
−
11
5.170
×
10
−
7
5.5
×
10
−
3
50
5.733
×
10
−
17
1.730
×
10
−
15
5.044
×
10
−
12
4.156
×
10
−
8A comparison of the HAM exact solutions and the analytical solutions for the electric potential field is shown in Figure3(a). Note that the solutions obtained through both techniques are smooth throughout the domain. Furthermore, the exact solution calculated through the HAM technique is in complete agreement with the analytical solution from equation (42). Similarly, Figure 3(b) presents a comparison of the HAM velocity solution with the analytical solutions obtained by Sun et al. [55], which is given as(43)Wη=−η2Γ4+2E¯sG1κ2I0κη−1I0κ+C1logη+C2.Figure 3
Comparison of the HAM results with the analytical approach for the solutions of electric potential profileψη[33] and velocity profile Wη [55] against η when R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G2=G3=G4=G5=1. Line: HAM exact solution. Delta: analytical solution.
(a)
(b)After using the boundary conditions (13) and (14), equation (43) becomes(44)Wη=−η2Γ4+2E¯sG1κ2I0κη−1I0κ+Γ4+2E¯sG11−I0κκ2I0κ.The values for unknown constantsΓ and E¯s obtained through the HAM formulation are used with κ=G1=1. It can be seen that the 50th-order HAM solution is fully consistent with the analytical solution.


### 5.2. Analysis and Discussion

The electroosmotic parameter κ depends on the pipe radius and Debye–Hückel parameter. An increase in this parameter implies an increase in the channel radius and a decrease in the width of EDL compact layer. Because of this change, the electrostatic potential in the pipe is expected to decrease. Additionally, EDL has minimal effects in channels with larger radius. Therefore, the electrostatic potential is expected to be immune to a change inκ when this parameter is significantly large. It can be seen in Figure 3(a) that the electrostatic potential decreases in the mid part of the pipe as κ increases. It can also be seen that the effect on electrostatic potential decreases as κ increases and it becomes negligible when κ in increased to a sufficiently large value. This trend is consistent with the study conducted by Mala et al. [33], which analyzed the Boltzmann–Poisson model for EDL.The two unknown constants, streaming potentialE¯s and pressure gradient Γ, are computed during the HAM computations. Figure 4 illustrates the trend in these quantities with respect to the electroosmotic parameter κ. Gong and Wu [56] report the streaming potential as the electric resistance number and state that the larger value of this constant means higher resistance due to EDL. Note that, as κ increases, the streaming potential E¯s decreases; but, this decrease becomes negligible after a significantly large value of κ as the E¯s approaches to zero. The effects of this decrease can be seen in the velocity and temperature profiles hereinafter. It shows that, as the EDL thickness decreases, the EDL effects decrease and become negligible eventually. In Figure 4(b), it can be seen that the pressure gradient constant starts decreasing rapidly from 12; but, as κ increases to a significantly large value, this change becomes negligible and eventually becomes constant. Physically, κ is inversely proportional to the thickness of EDL compact layer; so, after a significantly large value of κ, the radius does not change enough to influence the pressure gradient constant.Figure 4
The unknown constants: the streaming potentialE¯s and the pressure gradient constant Γ against the electroosmotic parameter κ when R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)Malvandi and Ganj [57] noted that the nanoparticle concentration in the mid part of the channel decreased with increase in thermophoretic parameter. They drew a similar conclusion in case of microtubes in another study [59]. The particle acceleration is higher for higher thermophoretic diffusion which causes a drop in the fluid nanoparticle concentration. Figure 5 illustrates the change in the nanoparticle volume fraction s(η) when the thermophoretic diffusion parameter and the Brinkman number are varied. Note in Figure 5(a) that increasing the thermophoretic diffusion in the flow causes a decrease in the nanoparticle concentration throughout the domain. Figure 5(b) shows change in the nanoparticle concentration when the Brinkman number is varied. Note that the nanoparticle concentration decreases with increase in Br due to the increased fluid viscosity.Figure 5
The nanoparticle profilesη against η with variation in the thermophoresis parameter (Nt) and the Brinkman number (Br) when R=1, ζ=1, κ=1, Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)The electrostatic potential is expected to affect the flow velocity in such a way that the resistance due to EDL will initially increase near the pipe boundary. However, this resistance will start to decrease asκ exceeds the max effect value κ∗. This behavior can be observed in Figure 6. The dimensionless skin friction coefficient is given as Cf=2W′1/Re, where Re=RρWm/μ is the Reynolds number. It can be seen in Figure 6(a) that the resistance on the pipe boundary initially increases with increase in κ, but this trends changes as κ exceeds κ∗=2.30675. κ∗ is obtained through maximization of Cf with respect to κ.Figure 6
The local skin friction coefficientCf and velocity profile Wη with variation in the electroosmotic parameter κ when R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)The velocityWη is observed against η with variation in κ in Figure 6(b). Note that there is a small change in the velocity profile when κ is changed. This change is maximum at κ=κ∗ and smaller for larger values of κ, which indicates the decrease of EDL effects with increase in the pipe radius; the velocity near the boundary decreases with increase in κ, but this decrease gets less prominent for higher values of κ and eventually, the velocity profile approaches to the case with no EDL. This is consistent with the results obtained by Gong and Wu [56], according to which the resistance due to EDL is greater near the boundaries as compared to the mid part of the channel. Additionally, according to Zhao et al. [54], the resistance in the velocity due to EDL is reduced with increase in κ, which is also evident in this case.Similarly, the electrostatic potential profile will affect the temperature distribution due to Joule heating effect in the energy equation as illustrated in Figure7. The dimensionless Nusselt number, given as Nu=−θ1, determines the dominance of convective heat transfer through the pipe/channel boundary. Plotted in Figure 7(a), the Nusselt number increases with increase in κ. However, this increase sufficiently slows down as κ>κ∗ and eventually becomes negligible when κ is sufficiently large. The temperature profile θη with change in κ is observed in Figure 7(b). Note that the temperature throughout the pipe decreases as κ increases. Since the EDL effects decrease with increase in κ, the effects from this term also decrease as a consequence, and this change becomes negligibly small as κ increases to a sufficiently large value.Figure 7
The local Nusselt number and velocity profileWη against η with variation in the electroosmotic parameter κ when R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)It is important to analyze the influence of the thermophoretic diffusion, Brownian diffusion, and viscosity on the temperature profile in a microradius flow. An increase in the thermophoretic diffusion in the flow can cause an increase in the fluid-particle acceleration, which as a consequence can cause an increase in the flow temperature profile. This trend was observed by Malvandi and Ganj [57] in their study of water/alumina nanofluid in a vertical microchannel. Their results state that when the thermophoretic parameter is increased, the temperature in the channel is increased. Similarly, there is a higher amount of molecular collision, which causes rise in the flow temperature, when Brownian diffusion in a flow is increased. This was observed by Animasaun et al. [58] in detail. According to their hypothesis, the flow temperature is bound to increase with an increase in the Brownian motion in the flow. The Brinkman number (Br) is the ratio of the heat generation due to viscosity to external heating. So, an increase in the Brinkman number implies an increase in the fluid viscosity, which can raise the temperature in the fluid. It shall be noted in Figure 8(a) that, as the thermophoretic parameter Nt increases, the temperature profile increases throughout the domain. Similarly, the temperature distribution enhances with an increase in the Brownian diffusion parameter, which is shown in Figure 8(b). Similarly, the temperature increases with increase in the Brinkman number, as shown in Figure 8(c). All these results are consistent with the previous studies on heat convection in fluid flow.Figure 8
The temperature profileθη against η with variation in the thermophoresis parameter (Nt), Brownian diffusion parameter (Nb), and the Brinkman number (Br) when R=1, ζ=1, κ=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)
(c)Figures9 and 10 present an entropy analysis for the flow. Note in equation (20) that the flow entropy is highly dependent on the temperature, velocity, and nanoparticle concentration of the flow. From the former analysis, it is obvious that the electroosmotic parameter κ effects the flow velocity and temperature fields. The effects on the entropy due to change in κ are illustrated in Figure 9. Note that, in Figure 9(a), the total entropy of the flow is maximum when the value of κ is near zero. The total flow entropy decreases with increase in κ. From previous analysis, it is evident that the EDL effect is maximum at κ=κ∗. This effect starts to decrease as κ>κ∗, which is evident in the flow entropy trend shown in the figure. The observed decrease is due to the decrease in temperature profile with increase in κ, which can be further understood by observing Figure 9(b), where the entropy field of the flow is observed at different values of κ. It can be seen that, for lower values of κ, i.e., κ<κ∗, there is more disturbance in the entropy near the pipe wall and the entropy is relatively higher. It reduces throughout the channel when κ is increased, but after a sufficiently large value of κ, this change becomes negligible.Figure 9
Change in the flow entropy with evolvingκ when R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)Figure 10
Total entropy of the flow with change in the Brinkman number (Br) and the thermophoretic diffusion parameter (Nt) when R=1, ζ=1, κ=1, Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.The effects of the thermophoretic diffusion parameter (Nt) and the Brinkman number (Br) on the total entropy of the flow can be seen in Figure 10. The total entropy of the channel increases with increase in Br, whereas this trend is reversed in case of Nt where it decreases with increase in the thermophoretic parameter. The reason for this trend is illustrated in Table 2. It shall be noted that, with increase in Br, there is a notable increase in the entropy due to heat transfer and viscous friction, whereas a minute decrease can be observed due to the nanoparticle and Joule heating. Hence, the total entropy of the flow is increased. On the other hand, entropy due to nanoparticle has a prominent decrease when Nt is increased; this causes a decrease in the total entropy of the flow. This can also be seen in Figures 5(a) and 8(a) where the decrease in the nanoparticle concentration (due to Nt) is greater compared to the increase in the temperature field.Table 2
Change in the total entropy and entropy terms due to the Brinkman number and the thermophoretic parameter.
Br
N
t
0.1
0.3
0.6
0.9
0.05
0.1
0.150
0.2
S
G
1.768
2.175
2.777
3.368
1.831
1.768
1.709
1.656
S
H
0.535
0.689
0.953
1.253
0.529
0.535
0.542
0.548
S
C
0.540
0.517
0.473
0.416
0.609
0.540
0.477
0.419
S
J
0.534
0.509
0.478
0.450
0.533
0.532
0.531
0.530
S
F
0.158
0.459
0.871
1.246
0.159
0.158
0.158
0.158


## 5.1. Solution Validation

The convergence of the HAM technique has been tested and is presented in Figure2. The 20th-order HAM solutions are used to determine the convergence domains for the parameters. It shall be noted that, for every convergence control parameter, there exists a convergence domain. Therefore, it can be stated that the applied technique is convergent. For the validation of the obtained solution, the approach used by [54] is used. Accuracy of the solution is checked using the maximum total average squared error which is given as(40)Em=maxEΨm,EWm,Eθm,Esm,where m represents the order of the HAM computation and EΨm,EWm,Eθm,andEsm denote the functions for error evaluation, which are defined, based on equations (9)–(12), given as(41)EΨm=∫01Ψ″+Ψ′η−κ2Ψ2dη,EWm=∫01W″+W′η+Γ−2G1E¯sΨ2dη,Eθm=∫01θ″+θ′η+Ntθ′2+Nbs′θ′+BrW′2+G2E¯s2+A1W+G3NbA2+G4NtA1,Esm=∫01s″+s′η+NtNbθ″+θ′η+A2W2dη.Figure 2
HAM convergence control parameter curves. Here,R=1, ζ=1, κ=1, Nb=Nt=1, Br=1, A1=A2=1, and G1=G2=G3=G4=G5=1.Equation (40) is used to compute the mth-order HAM computation error. Table 1 illustrates the errors computed using the following parameter values: R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1. Note in the table that the error decreases as the HAM computational order is increased; this shows that the results are convergent. It should be noted that the homotopy-Padé technique [50] was used to increase the convergence rate and obtain more accuracy, indicating the insignificance of the convergence control parameters which can hence be given any nonzero value. The 50th-order HAM solutions have a good enough accuracy to be used for the foregoing analysis. The HAM approach for solution computation can be further verified using the analytical solution obtained for equation (9). After application of the boundary conditions (15) and (16) and neglecting the imaginary parts, the solution is obtained as(42)ψη=ζI0κrI0κ,where I0 is the zeroth-order modified Bessel function of the first kind.Table 1
The maximum errorEm.
Order
κ
=
1
κ
=
3
κ
=
5
κ
=
10
10
41.766
8.94
×
10
6
2.08
×
10
6
794.505
20
0.52
0.025
8.063
3.651
30
5.138
×
10
−
7
4.803
×
10
−
7
8.2
×
10
−
4
2.557
40
5.877
×
10
−
12
5.552
×
10
−
11
5.170
×
10
−
7
5.5
×
10
−
3
50
5.733
×
10
−
17
1.730
×
10
−
15
5.044
×
10
−
12
4.156
×
10
−
8A comparison of the HAM exact solutions and the analytical solutions for the electric potential field is shown in Figure3(a). Note that the solutions obtained through both techniques are smooth throughout the domain. Furthermore, the exact solution calculated through the HAM technique is in complete agreement with the analytical solution from equation (42). Similarly, Figure 3(b) presents a comparison of the HAM velocity solution with the analytical solutions obtained by Sun et al. [55], which is given as(43)Wη=−η2Γ4+2E¯sG1κ2I0κη−1I0κ+C1logη+C2.Figure 3
Comparison of the HAM results with the analytical approach for the solutions of electric potential profileψη[33] and velocity profile Wη [55] against η when R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G2=G3=G4=G5=1. Line: HAM exact solution. Delta: analytical solution.
(a)
(b)After using the boundary conditions (13) and (14), equation (43) becomes(44)Wη=−η2Γ4+2E¯sG1κ2I0κη−1I0κ+Γ4+2E¯sG11−I0κκ2I0κ.The values for unknown constantsΓ and E¯s obtained through the HAM formulation are used with κ=G1=1. It can be seen that the 50th-order HAM solution is fully consistent with the analytical solution.


## 5.2. Analysis and Discussion

The electroosmotic parameter κ depends on the pipe radius and Debye–Hückel parameter. An increase in this parameter implies an increase in the channel radius and a decrease in the width of EDL compact layer. Because of this change, the electrostatic potential in the pipe is expected to decrease. Additionally, EDL has minimal effects in channels with larger radius. Therefore, the electrostatic potential is expected to be immune to a change inκ when this parameter is significantly large. It can be seen in Figure 3(a) that the electrostatic potential decreases in the mid part of the pipe as κ increases. It can also be seen that the effect on electrostatic potential decreases as κ increases and it becomes negligible when κ in increased to a sufficiently large value. This trend is consistent with the study conducted by Mala et al. [33], which analyzed the Boltzmann–Poisson model for EDL.The two unknown constants, streaming potentialE¯s and pressure gradient Γ, are computed during the HAM computations. Figure 4 illustrates the trend in these quantities with respect to the electroosmotic parameter κ. Gong and Wu [56] report the streaming potential as the electric resistance number and state that the larger value of this constant means higher resistance due to EDL. Note that, as κ increases, the streaming potential E¯s decreases; but, this decrease becomes negligible after a significantly large value of κ as the E¯s approaches to zero. The effects of this decrease can be seen in the velocity and temperature profiles hereinafter. It shows that, as the EDL thickness decreases, the EDL effects decrease and become negligible eventually. In Figure 4(b), it can be seen that the pressure gradient constant starts decreasing rapidly from 12; but, as κ increases to a significantly large value, this change becomes negligible and eventually becomes constant. Physically, κ is inversely proportional to the thickness of EDL compact layer; so, after a significantly large value of κ, the radius does not change enough to influence the pressure gradient constant.Figure 4
The unknown constants: the streaming potentialE¯s and the pressure gradient constant Γ against the electroosmotic parameter κ when R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)Malvandi and Ganj [57] noted that the nanoparticle concentration in the mid part of the channel decreased with increase in thermophoretic parameter. They drew a similar conclusion in case of microtubes in another study [59]. The particle acceleration is higher for higher thermophoretic diffusion which causes a drop in the fluid nanoparticle concentration. Figure 5 illustrates the change in the nanoparticle volume fraction s(η) when the thermophoretic diffusion parameter and the Brinkman number are varied. Note in Figure 5(a) that increasing the thermophoretic diffusion in the flow causes a decrease in the nanoparticle concentration throughout the domain. Figure 5(b) shows change in the nanoparticle concentration when the Brinkman number is varied. Note that the nanoparticle concentration decreases with increase in Br due to the increased fluid viscosity.Figure 5
The nanoparticle profilesη against η with variation in the thermophoresis parameter (Nt) and the Brinkman number (Br) when R=1, ζ=1, κ=1, Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)The electrostatic potential is expected to affect the flow velocity in such a way that the resistance due to EDL will initially increase near the pipe boundary. However, this resistance will start to decrease asκ exceeds the max effect value κ∗. This behavior can be observed in Figure 6. The dimensionless skin friction coefficient is given as Cf=2W′1/Re, where Re=RρWm/μ is the Reynolds number. It can be seen in Figure 6(a) that the resistance on the pipe boundary initially increases with increase in κ, but this trends changes as κ exceeds κ∗=2.30675. κ∗ is obtained through maximization of Cf with respect to κ.Figure 6
The local skin friction coefficientCf and velocity profile Wη with variation in the electroosmotic parameter κ when R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)The velocityWη is observed against η with variation in κ in Figure 6(b). Note that there is a small change in the velocity profile when κ is changed. This change is maximum at κ=κ∗ and smaller for larger values of κ, which indicates the decrease of EDL effects with increase in the pipe radius; the velocity near the boundary decreases with increase in κ, but this decrease gets less prominent for higher values of κ and eventually, the velocity profile approaches to the case with no EDL. This is consistent with the results obtained by Gong and Wu [56], according to which the resistance due to EDL is greater near the boundaries as compared to the mid part of the channel. Additionally, according to Zhao et al. [54], the resistance in the velocity due to EDL is reduced with increase in κ, which is also evident in this case.Similarly, the electrostatic potential profile will affect the temperature distribution due to Joule heating effect in the energy equation as illustrated in Figure7. The dimensionless Nusselt number, given as Nu=−θ1, determines the dominance of convective heat transfer through the pipe/channel boundary. Plotted in Figure 7(a), the Nusselt number increases with increase in κ. However, this increase sufficiently slows down as κ>κ∗ and eventually becomes negligible when κ is sufficiently large. The temperature profile θη with change in κ is observed in Figure 7(b). Note that the temperature throughout the pipe decreases as κ increases. Since the EDL effects decrease with increase in κ, the effects from this term also decrease as a consequence, and this change becomes negligibly small as κ increases to a sufficiently large value.Figure 7
The local Nusselt number and velocity profileWη against η with variation in the electroosmotic parameter κ when R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)It is important to analyze the influence of the thermophoretic diffusion, Brownian diffusion, and viscosity on the temperature profile in a microradius flow. An increase in the thermophoretic diffusion in the flow can cause an increase in the fluid-particle acceleration, which as a consequence can cause an increase in the flow temperature profile. This trend was observed by Malvandi and Ganj [57] in their study of water/alumina nanofluid in a vertical microchannel. Their results state that when the thermophoretic parameter is increased, the temperature in the channel is increased. Similarly, there is a higher amount of molecular collision, which causes rise in the flow temperature, when Brownian diffusion in a flow is increased. This was observed by Animasaun et al. [58] in detail. According to their hypothesis, the flow temperature is bound to increase with an increase in the Brownian motion in the flow. The Brinkman number (Br) is the ratio of the heat generation due to viscosity to external heating. So, an increase in the Brinkman number implies an increase in the fluid viscosity, which can raise the temperature in the fluid. It shall be noted in Figure 8(a) that, as the thermophoretic parameter Nt increases, the temperature profile increases throughout the domain. Similarly, the temperature distribution enhances with an increase in the Brownian diffusion parameter, which is shown in Figure 8(b). Similarly, the temperature increases with increase in the Brinkman number, as shown in Figure 8(c). All these results are consistent with the previous studies on heat convection in fluid flow.Figure 8
The temperature profileθη against η with variation in the thermophoresis parameter (Nt), Brownian diffusion parameter (Nb), and the Brinkman number (Br) when R=1, ζ=1, κ=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)
(c)Figures9 and 10 present an entropy analysis for the flow. Note in equation (20) that the flow entropy is highly dependent on the temperature, velocity, and nanoparticle concentration of the flow. From the former analysis, it is obvious that the electroosmotic parameter κ effects the flow velocity and temperature fields. The effects on the entropy due to change in κ are illustrated in Figure 9. Note that, in Figure 9(a), the total entropy of the flow is maximum when the value of κ is near zero. The total flow entropy decreases with increase in κ. From previous analysis, it is evident that the EDL effect is maximum at κ=κ∗. This effect starts to decrease as κ>κ∗, which is evident in the flow entropy trend shown in the figure. The observed decrease is due to the decrease in temperature profile with increase in κ, which can be further understood by observing Figure 9(b), where the entropy field of the flow is observed at different values of κ. It can be seen that, for lower values of κ, i.e., κ<κ∗, there is more disturbance in the entropy near the pipe wall and the entropy is relatively higher. It reduces throughout the channel when κ is increased, but after a sufficiently large value of κ, this change becomes negligible.Figure 9
Change in the flow entropy with evolvingκ when R=1, ζ=1, Br=1, Nt=Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.
(a)
(b)Figure 10
Total entropy of the flow with change in the Brinkman number (Br) and the thermophoretic diffusion parameter (Nt) when R=1, ζ=1, κ=1, Nb=1, A1=A2=1, and G1=G2=G3=G4=G5=1.The effects of the thermophoretic diffusion parameter (Nt) and the Brinkman number (Br) on the total entropy of the flow can be seen in Figure 10. The total entropy of the channel increases with increase in Br, whereas this trend is reversed in case of Nt where it decreases with increase in the thermophoretic parameter. The reason for this trend is illustrated in Table 2. It shall be noted that, with increase in Br, there is a notable increase in the entropy due to heat transfer and viscous friction, whereas a minute decrease can be observed due to the nanoparticle and Joule heating. Hence, the total entropy of the flow is increased. On the other hand, entropy due to nanoparticle has a prominent decrease when Nt is increased; this causes a decrease in the total entropy of the flow. This can also be seen in Figures 5(a) and 8(a) where the decrease in the nanoparticle concentration (due to Nt) is greater compared to the increase in the temperature field.Table 2
Change in the total entropy and entropy terms due to the Brinkman number and the thermophoretic parameter.
Br
N
t
0.1
0.3
0.6
0.9
0.05
0.1
0.150
0.2
S
G
1.768
2.175
2.777
3.368
1.831
1.768
1.709
1.656
S
H
0.535
0.689
0.953
1.253
0.529
0.535
0.542
0.548
S
C
0.540
0.517
0.473
0.416
0.609
0.540
0.477
0.419
S
J
0.534
0.509
0.478
0.450
0.533
0.532
0.531
0.530
S
F
0.158
0.459
0.871
1.246
0.159
0.158
0.158
0.158


## 6. Conclusions

A fully developed nanofluid flow in a microradius pipe with the presence of electrical effects from an electrical double layer has been studied. The governing equations embodying, the Poisson–Boltzmann equation, momentum equation, energy equation, and Buongiorno’s nanofluid model equation, have been transformed using the appropriate similarity transformations. The unknown pressure gradient constant has been computed using the average velocity condition, and the streaming potential due to EDL has derived for the considered problem. Highly accurate exact solutions have been obtained using the homotopy analysis method. The findings of this study are as follows.The model used for formulation of the nanofluid flow is consistent with the experimental studies in the past. Change in the electroosmotic parameter has high influence on the flow when it is smaller than the max effect valueκ∗, which was found to be 2.30675. The EDL effects reduce and eventually become negligible as κ>κ∗. The electrostatic potential, velocity, temperature, and nanoparticle volume fraction are affected by the evolution in viscosity, thermophoretic diffusion, Brownian diffusion, and electroosmosis. The flow entropy increases with increase in viscosity and thermophoretic diffusion, while it decreases when the EDL compact layer is reduced. Furthermore, the nanoparticle and temperature concentration are found to be two of the main factors causing the change in flow entropy.



---

*Source: 1723256-2020-02-18.xml*