Thermal and Flow Analysis of Fully Developed Electroosmotic Flow in Parallel-Plate Micro- and Nanochannels with Surface Charge-Dependent Slip

Abstract

This study analytically investigates the coupled effects of surface charge and boundary slip on the fully developed electroosmotic flow and thermal transfer in parallel plate micro and nanochannels under the high zeta potential. The electric potential, velocity, temperature, flow rate, and Nusselt number are obtained analytically. The main results are that the velocity of bulk flow is significantly reduced in the presence of the surface charge-dependent slip. Moreover, the maximum velocity at ζ = −125 mV is approximately twice as large as that at ζ = −25 mV. The velocity and dimensionless temperature increase as the zeta potential increases. The dimensionless temperature of the surface charge-dependent slip flow is larger than that of the surface charge-independent slip flow. For the surface charge-dependent slip flow, the maximum temperature at ζ = −125 mV is approximately four times larger than that at ζ = −25 mV. The Nusselt number decreases with Joule heating and increases with a positive heat transfer coefficient. The Nusselt number decreases as the electric field and the magnitude of the zeta potential increase. In the surface charge-dependent slip flows, the Nusselt number is smaller than that in the surface charge-independent slip flows.

1. Introduction

Over the last decades, the micro- and nanoscale fluid systems based on physics, medicine, chemistry and biology have been widely applied in various fields [1,2], such as the separation and detection of chemical and biological samples, the design of heat and mass transfer system, micro-and nanoflow control in lab-on-chip [2]. In these applications, the thermal properties of fluid flows are significant factors affecting the performance of micro-nanoscale fluid systems, including the temperature of fluid could affect the properties of biological and medical samples [3], as well as the formation and stability of nanobubbles at the solid–liquid interface [4]. Traditionally, compared with macroscopic flows, the reduced characteristic scale of microfluidic channels results in surface forces that are much larger than the volume force. Therefore, for certain aspects like surface force, capillary effect, slip effect, and rapid heat transfer effect [5], the micro- and nanofluidic mechanics exhibit new characteristics distinguished from the macro-scale fluid flow. At present, several mechanisms for fluid-driving in micro- and nanoflows have been put forward in the existing literature, such as the micropumps by pressure gradient [6], external electric field, external magnetic field [7,8,9,10], surface tension, and acoustic wave [11]. In most cases, a solid surface acquires a negative electric charge when brought into contact with an electrolyte solution. In return, the charged solid surface redistributes nearby ions in the surrounding solution forming a double-layer structure of ions in the vicinity of the surface of the solid, which is referred to as the electric double layer (EDL) [2,12]. The EDL is a region close to the charged surface, where the counter ions exceed the ions, thus neutralizing the surface charge. When an axial electric field is applied, the electric power generated by the interaction between the net charge density in EDL and the applied electric field can drive the flow. This kind of flow is known as electroosmotic flow (EOF) [2,12].

Due to its relative advantages, such as simple design requirements and no need for moving elements, the electroosmosis force has been extensively discussed as a more effective driving mechanism than a simple pressure driving mechanism. In this respect, much of the literature on the microchannel electroosmosis flow by Newtonian fluid and non-Newtonian fluid has been well studied [13,14,15,16,17,18,19,20,21,22,23,24,25]. Yang et al. [17] analyzed the EOF in an annulus with high zeta potentials. Kadet and Koryuzlov [22] have presented the results of experimental studies of EOF of two immiscible liquids in a thin slit channel. Through experimental, simulation and theoretical analyses in polymer-coated slits, Monteferrante and Melchionna [23] pointed out that electroosmotic suppression arises from the frictional forces acting on the ionic currents.

However, when the voltage is high, Joule heating is inevitable in the process of EOF. How to dramatically improve the performance of conventional heat transfer has become a major challenge for scientists and engineers. Therefore, heat transfer has become the focus of EOF research in recent years [26,27,28,29,30,31,32,33,34,35,36,37,38]. Under the condition of a constant heat flux boundary, Maynes and Webb [31] studied the thermal transport of EOF in parallel plates and circular microchannels and obtained the analytical expression of temperature distribution. In addition, they also studied the viscous dissipation effect [32]. Liechty et al. [33] investigated the fully-developed convection heat transfer of EOF in a circular microtube with arbitrary wall zeta potentials under imposed wall temperatures and imposed wall heat fluxes. Tang et al. [34] have obtained a numerically transient temperature field that takes into account Joule heating in a microchannel. Chakraborty [35] has gained analytical solutions of the temperature and the Nusselt number of the fully developed flow driven by both electroosmotic force and imposed pressure gradient. Sadeghi and Saidi [36] have investigated the influence of viscous dissipation on thermal transport characteristics of the fully developed combined electroosmotic and pressure-driven flow in a parallel-plate microchannel with uniform wall heat flux. The research results show that viscous dissipation and Joule heating dominate the heat transfer problem in EOF.

However, most studies on EOF problems are limited to no-slip boundary conditions. Actually, the experimental results show that when the fluid flows over the solid surface, there may be relative motion between the solid surface and the liquid near the solid surface. The velocity of the solid surface is called slip velocity u_slip, which is related to slip length b. The slip length is usually in the range of several nanometers to tens of micrometers [39]. Previous studies have shown that EDL and slip can significantly affect micro- and nanoscale fluid behaviors [3,40,41,42,43,44,45,46,47,48,49]. In the case of electroosmotic flow, Rojas et al. [19] conducted a dispersion analysis in the neutral non-reacting solute in the circular microchannel with boundary slip and concluded that the existence of slip can greatly magnify the dispersion increase caused by the induced pressure gradient caused by the non-uniform wall potential [20]. Kundu and Saha [24] discussed the thermal transport behavior of electrostatic-driven flow in microchannels based on the no-slip, first- and second-order slip at the boundary of velocity distribution and the no-jump, first- and second-order heat slip of thermal response under the condition of maintaining uniform wall heat flux. In a microchannel between two parallel plates with imposed heat flux, the liquid flows with the slip boundary condition have been numerically investigated by Ngoma and Erchiqui [38]. Park and Kim [46] studied the Electrokinetic slip flow through hydrophobic microchannels and obtained an analytical solution of volumetric flow rate and solute retention time.

Although most of the studies in the literature were made to investigate the effects of EOF and convective heat transfer on liquid flow, the impacts of the surface charge about boundary slip on heat transfer in electroosmosis flows were rarely concerned about [40,50]. Under the influence of the surface charge, EDL with an opposite net charge is formed in the solution near the surface of the wall. Then, this charge distribution near the wall results in an attractive electric force between the charged wall and the nearby fluid. Eventually, this interaction between the wall and the fluid further influences the slip boundary conditions [41,50]. Utilizing the molecular dynamics simulation method, Joly et al. [50] studied and used a mathematical model to describe the effect of surface charge on the boundary slip. Jing and Bhushan [41,51] used atomic force microscopy to study the relationship between the boundary slip and the surface charge and found the same relationship as that obtained by Joly et al. [50]. So far, although most researchers have found the dependence of boundary slip on a surface charge, the study of fluid types is mainly limited to Newtonian fluid [50,51,52,53,54,55,56]. Bhushan et al. [40,41] studied the effect of surface charge-dependent boundary slip on the electroviscous effect between microchannels. Pan et al. [57,58] further investigated the effect of surface charge-dependent slip on the electroviscous effect and heat transfer in a microchannel and obtained analytical expressions for velocity, temperature distributions, and Nusslt number.

The EOFs in slit channels have a variety of practical applications, such as micro/nano-pumping, mixing, and oil extraction [22]. However, a large number of previous studies rarely discussed this problem theoretically. The purpose of this paper is to investigate the effects of the surface charge-dependent slip on the combined electroosmotic and pressure-driven flow in a parallel-plate micro and nanochannel with high zeta potential and constant wall heat flux. The paper is organized as follows: Section 2 presents the mathematical modeling. Results and comparisons are discussed in Section 3. The conclusions are drawn in Section 4.

2. Mathematical Modeling

2.1. Electric Potential Distribution

In this study, the steady combined electroosmotic and pressure-driven flow of a symmetric binary electrolyte between two parallel plates with surface charge-dependent slip and high zeta potential is studied with density ρ, dynamic viscosity μ, and electric conductivity σ of the fluid. A Cartesian orthonormal coordinate system (x, y, z) is used with the origin fixed at the middle of the channel, as shown in Figure 1. The channel width W in the z-direction and length L in the x-direction is much larger than the channel height 2H in the y-direction, i.e., 2H << L and W << L. EOF in the x-direction is generated when an electric field is applied along the x-direction. On the boundary, the effect of surface charge density σ_s on the boundary slip b is considered [41,53]:

$$b=\frac{b_0}{1+\left(1/\alpha \right){\left({\sigma }_s{d}^2 /e\right)}^2\left(l_B/d^2\right)b_0},$$

(1)

where b₀ is the original slip length without considering surface charge effect, α~1 is a numerical factor, e is the elementary charge, d is the equilibrium distance of Lennard-Jones potential, and l_B = e²/(4πεk_BT_av), ε is the permittivity of the electrolyte solution, k_B is the Boltzmann constant and T_av is the average absolute temperature over the channel cross-section [53,54]. According to the EDL theory, the electric potential $\psi$ and the net charge density ρ_e satisfy Poisson equation [2,12]

$${\nabla }^2\psi =-\frac{{\rho }_e}{\epsilon }.$$

(2)

The ideal solution of fully dissociated symmetric salt, the electric charge density is given by [2,12]:

$${\rho }_e=-2n_0{z}_{v}e\sinh \frac{z_{v}e\psi }{k_B{T}_{av}},$$

(3)

where n₀ is the ion density, z_v is the valence number of ions in solution. Using Equations (2) and (3), we can obtain the Poisson–Boltzmann equation

$$\frac{{\mathrm{d}}^2\psi }{\mathrm{d}{y}^2}=\frac{2n_0{z}_{v}e}{\epsilon }\sinh \left(\frac{z_{v}e\psi }{k_B{T}_{av}}\right) .$$

(4)

The following dimensionless parameters are introduced:

$$\overline{\psi }=\frac{z_{v}e\psi }{k_B{T}_{av}} , \overline{y}=\frac{y}{H} , K=\kappa H,$$

(5)

where K is called the nondimensional electrokinetic width, κ = (2n₀z_v²e²/(εk_BT_av))^1/2 is Debye–Hückel parameter and 1/κ is the characteristic thickness of the EDL. The dimensionless form of the Poisson–Boltzmann equation can be expressed as:

$$\frac{{\mathrm{d}}^2\overline{\psi }}{\mathrm{d}{\overline{y}}^2}=K^2\sinh \overline{\psi }.$$

(6)

Equation (6) is subject to the following boundary conditions:

$${\overline{\psi }|}_{\overline{y}=-1}=Z, \frac{\mathrm{d}\overline{\psi }}{\mathrm{d}\overline{y}}{|}_{\overline{y}=0}=0,$$

(7)

here, Z is the dimensionless zeta potential, i.e., Z = z_veζ/(k_BT_av).

From Equation (6) and boundary conditions (7), one obtains the dimensionless potential distribution:

$$\overline{\psi }=4{\tanh }^{-1}\left[e^{-K\left(1+\overline{y}\right)}\tanh \left(Z/4\right)\right].$$

(8)

2.2. Velocity Distribution

Since the velocity of EOF and the hydraulic diameter of micro/nano channels are very small, the Reynolds number is very low, and hence fully-developed laminar flow can be assumed. In the fully-developed parallel flow $(\frac{\partial u}{\partial x}=0,v=0,w=0$), the Navier–Stokes equation becomes linear. Here, u, v, and w are the velocities in the x, y, and z-directions, respectively. The momentum equation in the x-direction under the combined action of pressure and electric field E₀ is as follows:

$$\mu \frac{{\mathrm{d}}^{2}u}{\mathrm{d}{y}^2}=\frac{\mathrm{d}p}{\mathrm{d}x}-{\rho }_e{E}_0.$$

(9)

Equation (9) is subject to the following boundary conditions:

$$u{|}_{y=-H}=b\frac{\mathrm{d}u}{\mathrm{d}y}{|}_{y=-H}, \frac{\mathrm{d}u}{\mathrm{d}y}{|}_{y=0}=0,$$

(10)

where the slip length b is given by Equation (1).

According to the charge conservation law, the surface charge should be equal to the net charge in the fluid. Therefore

$${\sigma }_s=-{{\displaystyle \int }}_{-H}^0{\rho }_e\mathrm{d}y=-\epsilon \frac{\mathrm{d}\psi }{\mathrm{d}y}{|}_{y=-H}.$$

(11)

To write the governing equations and the boundary conditions in dimensionless form, let:

$$\overline{u}=\frac{u}{U_{eo}}, U_{eo}=-\frac{\epsilon \zeta E_0}{\mu },\overline{ b}=\frac{b}{H},$$

(12)

where U_eo is the steady Helmholtz–Smoluchowski velocity. First, substituting Equation (12) in the governing equations and in their boundary conditions and then solving the corresponding dimensionless boundary value problem by integrating Equation (9) twice, we obtain

$$\overline{u}=u_r\left(1+2\overline{b}-{\overline{y}}^2\right)+\left[Z+2\overline{b}K\sinh \left(Z/2\right)-\overline{\psi }\right]/Z,$$

(13)

in which u_r = u_max/U_eo is the ratio of the velocity of pressure-driven flow and that of EOF between two parallel plates, u_max = −H² (2μ)⁻¹ dp/dx. Note that u_r = 0 corresponds to a pure EOF.

According to the definition of the dimensional volumetric flow rate, we have

$$Q^{*}=2{{\displaystyle \int }}_0^W\mathrm{d}z{{\displaystyle \int }}_{-H}^{0}u\mathrm{d}y==2W{U}_{eo}H{{\displaystyle \int }}_{-1}^0\overline{u}\mathrm{d}\overline{y}.$$

(14)

Introduce a dimensionless volumetric flow rate

$$Q=\frac{Q^{*}}{W{U}_{eo}H}.$$

(15)

Substituting Equation (13) into Equation (14), the dimensionless volumetric flow rate can be expressed as

$$Q=4u_r\left(\frac{1}{3}+\overline{b}\right)+2+\frac{4\overline{b}K}{Z}\sinh (\frac{Z}{2})-\frac{2}{Z}{{\displaystyle \int }}_{-1}^0\overline{\psi }\mathrm{d}\overline{y}.$$

(16)

The integral in Equation (16) can be calculated numerically using Gaussian quadrature.

2.3. Temperature Distribution

Energy conservation, including the effects of Joule heating and viscous dissipation, can be written as

$$\rho c_{p}u\frac{\partial T}{\partial x}=k\left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}\right)+s+\mu {\left(\frac{\mathrm{d}u}{dy}\right)}^2,$$

(17)

where s and μ(du/dy)² denote the rate of volumetric heat generation due to Joule heating and viscous dissipation, respectively. s = i_e²σ_ρ, where i_e is the current density and σ_ρ is the liquid electrical resistivity. The current density is given by [33,36]

$$i_e=\frac{E_0}{{\sigma }_{\rho }}\cosh \overline{\psi },$$

(18)

where σ_ρ⁻¹ = 2z_v²e²Dn₀/(k_BT_av), D is the diffusivity of ions in the electrolyte. The boundary conditions for the energy equation are as follows:

$$\frac{\partial T}{\partial y}{|}_{y=0}=0,k\frac{\partial T}{\partial y}{|}_{y=-H}=-q or T{|}_{y=-H}=T_w\left(\mathrm{x}\right).$$

(19)

The dimensionless temperature θ is introduced

$$\theta \left(y\right)=\frac{T-T_w}{qH/k}.$$

(20)

This depends only on y for the fully developed flow. Differentiating both sides of Equation (20) with respect to x gives

$$\frac{\partial T}{\partial x}=\frac{\mathrm{d}{T}_w}{\mathrm{d}x}=\frac{\mathrm{d}{T}_b}{\mathrm{d}x}=const,$$

(21)

in which T_b is the bulk temperature. Considering an energy balance on a control volume with thickness dx along the centerline of the channel, we can write Equation (17) as:

$$\rho c_p\frac{\partial T_b}{\partial x}{{\displaystyle \int }}_{-H}^{0}u\mathrm{d}y=q+\frac{E_0^2}{{\sigma }_{\rho }}{{\displaystyle \int }}_{-H}^0{\cosh }^2\overline{\psi }\mathrm{d}y+\mu {{\displaystyle \int }}_{-H}^0{\left(\frac{\mathrm{d}u}{dy}\right)}^2\mathrm{d}y.$$

(22)

From Equation (22), we have

$$\frac{\mathrm{d}{T}_b}{\mathrm{d}x}=\frac{1}{\rho c_p{u}_{m}H}\left(q+\frac{H{E}_0^2}{{\sigma }_{\rho }}\beta +\frac{\mu U_{eo}^2}{H}\gamma \right),$$

(23)

where u_m, β, γ are given by

$$u_m=U_{eo}{{\displaystyle \int }}_{-1}^0\overline{u}\mathrm{d}\overline{y},\beta ={{\displaystyle \int }}_{-1}^0{\cosh }^2\overline{\psi }\mathrm{d}\overline{y}, \gamma ={{\displaystyle \int }}_{-1}^0{\left(\frac{\mathrm{d}\overline{u}}{\mathrm{d}\overline{y}}\right)}^2\mathrm{d}\overline{y}.$$

(24)

So the energy equation in dimensionless form can be written as

$$\frac{{\mathrm{d}}^2\theta }{\mathrm{d}{\overline{y}}^2}=\frac{1+S\beta +Br\gamma }{{\overline{u}}_m}\overline{u}-S{\cosh }^2\overline{\psi }-4Br{\left(u_r\overline{y}-\frac{K}{Z}\sinh \frac{\overline{\psi }}{2}\right)}^2,$$

(25)

where ${\overline{u}}_m$ = u_m/U_eo is the dimensionless axial mean velocity, Br = μU²_eo/qH is the Brinkman number, and S = E₀²H/σ_ρq is the dimensionless volumetric heat generation due to Joule heating. The thermal boundary conditions in the dimensionless form can be written as

$$\frac{\mathrm{d}\theta }{\mathrm{d}\overline{y}}{|}_{\overline{y}=0}=0, \theta {|}_{\overline{y}=-1}=0.$$

(26)

Integrating Equation (25) twice and using the boundary conditions (26), we obtain

$$\theta =f\left(\overline{y}\right)-f\left(-1\right),$$

(27)

here

$$\begin{array}{cc}f\left(\overline{y}\right)=\hfill & \frac{1+S\beta +Br\gamma }{2{\overline{u}}_m}\{u_r\left[\left(1+2\overline{b}\right){\overline{y}}^2-\frac{{\overline{y}}^4}{6}\right]-\frac{2}{Z}{{\displaystyle \int }}_0^{\overline{y}}\mathrm{d}\tilde{y}{{\displaystyle \int }}_0^{\tilde{y}}\widehat{\psi }\mathrm{d}\widehat{y}+\left[Z+2\overline{b}K\sinh \left(\frac{Z}{2}\right)\right]\hfill \\ & -S{{\displaystyle \int }}_0^{\overline{y}}\mathrm{d}\tilde{y}{{\displaystyle \int }}_0^{\tilde{y}}{\cosh }^2\widehat{\psi }\mathrm{d}\widehat{y}-4Br{{\displaystyle \int }}_0^{\overline{y}}\mathrm{d}\tilde{y}{{\displaystyle \int }}_0^{\tilde{y}}{\left[u_r\widehat{y}-\frac{K}{Z}\sinh \frac{\widehat{\psi }}{2}\right]}^2\mathrm{d}\widehat{y}.\hfill \end{array}$$

(28)

In the above equation, for simplicity, $\widehat{\psi }=\overline{\psi }\left(\widehat{y}\right)$. To obtain the Nusselt number, the dimensionless bulk temperature θ_b must be calculated first, which is given by

$${\theta }_b=\frac{{{\displaystyle \int }}_{-1}^0\overline{u}\theta \mathrm{d}\overline{y}}{{{\displaystyle \int }}_{-1}^0\overline{u} \mathrm{d}\overline{y}}=\frac{{{\displaystyle \int }}_{-1}^0\overline{u}\theta \mathrm{d}\overline{y}}{{\overline{u}}_m},$$

(29)

and Nusselt number is defined as

$$Nu=\frac{2hH}{k}=\frac{2qH}{k\left(T_w-T_b\right)}=-\frac{2}{{\theta }_b}.$$

(30)

3. Results and Discussion

In this section, the discrete numerical solutions of flow rate Q, dimensionless temperature θ and Nusselt number are obtained by Gaussian integration in two parallel plates with surface charge-dependent slip, and constant heat flux, together with the analytical and semi-analytic solutions for velocity, flow rate, and temperature. In addition, a series of results are displayed in the form of graphs to intuitively express the coupling effect of surface charge and slip on the combined pressure and electroosmotic-driven flow in a parallel-plate micro- and nanochannel. In order to obtain veritable and effective results, the following typical parametric values are used [2,12,53,54]: d = 0.4 × 10⁻⁹ m, α = 1, µ = 1.01 × 10⁻³ Pa∙s, ρ = 10³ kgm⁻³, e = 1.6 × 10⁻¹⁹ C, z = 1, T_av = 298 K, ε = 7 × 10⁻¹⁰ Fm⁻¹, k_B = 1.38 × 10⁻²³ JK⁻¹, D = 1.612 × 10⁻⁹ m²s⁻¹, without special instructions dp/dx = 0 Pa/m.

In Figure 2, the present results in the case of charge-independent slip are compared with those obtained by Mondal et al. [37] and Ngoma et al. [38]. Ngoma and Erchiqui [38] investigated the fluid flow with the slip boundary condition in a microchannel between two parallel plates; that is, the case of low zeta potential and the slip length b = b₀. In addition, Mondal and Shit [37] took into account the non-uniform walls varying sinusoidally, which can be simplified into the model in this paper in the limit case. The present results are found to be in agreement with those in Refs. [37,38]. Moreover, there is a slight error near the wall, which is mainly caused by the Debye–Hückel approximation used in Refs. [37,38].

In Figure 3, the non-slip velocity profile is plotted as a function of the distance for different values of the zeta potential and slip length. It can be seen from Figure 3a that for large zeta potentials, the velocity distribution of no-slip flow varies with zeta potential in the thin layer near the solid walls. As shown in Figure 3, it is found that the velocities of both the no-slip and slip flows increase with the zeta potential. For the surface charge-dependent slip flow, the maximum velocity at ζ = −125 mV is approximately twice as large as that at ζ = −25 mV. The study also demonstrated that the velocity of the surface charge-independent slip flow is larger than that of the surface charge-dependent flow. The surface charge-dependent expresses the case of b = b₀. In addition, the difference between the velocities of the surface charge-dependent and charge-independent slip flows increases with zeta potential, owing to the slip length decreasing with zeta potential, as shown in Equation (1).

In Figure 4, the effect of the surface charge-dependent slip on the flow rate is investigated. Here, the Q_ID represents the flow rate without considering the effect of the surface charges on the slip length. The surface charge effect is increased when the magnitude of zeta potential ζ increases and the slip length b₀ increases. Obviously, the slip length decreases with the zeta potential. A reasonable explanation is the increase of the attractive electrostatic force between the charged solid surface and the liquid near the solid–liquid interface. In addition, Figure 4 shows that the surface charge effect on the slip length in the microchannel is smaller than in the nanochannel.

In Figure 5, the temperature profile is plotted against the distance for different values of the zeta potential and the slip length. As seen in Figure 5, the dimensionless temperature in the nanochannel increases when the magnitude of the zeta potential increases. The reason is that as the magnitude of the zeta potential increases, the velocity and hence the effective Joule heating increases, and thus the heat transfer can be enhanced. Furthermore, it can be seen from Figure 5 that the dimensionless temperature in the surface charge-dependent slip flow is larger than that in the surface charge-independent slip flow. The viscous dissipation behaves like an energy source, increasing the temperature of the fluid, especially near the wall. The reason is that the highest shear rate occurs at this region while it is zero at the centerline. In addition, since the walls are cooling, heat is transferred from the interior to the exterior of the fluid, and fluid convection carries away a fraction of this heat, so the amplitude of the fluid velocity is large, further leading to a decrease in the amplitude of the temperature. For the surface charge-dependent slip flow, the maximum temperature at ζ = −125 mV is approximately four times larger than that at ζ = −25 mV.

The total Joule heating per unit length is defined as

$$E_{JH}={{\displaystyle \int }}_{-H}^H{i}_e^2{\sigma }_{\rho }\mathrm{d}y.$$

(31)

Nondimensionalizing this parameter by the total Joule heating per unit length in the case of low zeta potential, 2H (E₀²/σ_ρ), we obtain

$$E_e=\frac{E_{JH}}{2H{E}_0^2/{\sigma }_{\rho }}={{\displaystyle \int }}_{-1}^0{\cosh }^2\overline{\psi }\mathrm{d}\overline{y}=\beta .$$

(32)

The dimensionless total Joule heating per unit length of the channel E_e is plotted as a function of zeta potential ζ in Figure 6 for different values of H and n₀. It is found that for low zeta potential (ζ < 25 mV), E_e = 1 at all values of H. The reason can be illustrated that for low zeta potential, cosh $\overline{\psi }$~1. For a given zeta potential, the effective Joule heating E_e decreases as the half-height H of the micro and nanochannel increases. On the contrary, for a given H, the zeta potential strengths the effective Joule heating E_e, especially in nanochannels. By comparing Figure 6a,b, it can be found that ion density n₀ can inhibit the increase of Joule heat.

For a given H, as the zeta potential ζ increases, the effective Joule heating E_e increases, especially in nanochannels. In addition, as the ion density n₀ increases, the effective Joule heat E_e decreases. The reason is that an increase in the channel height H, or in the ion density n₀, implies an increase in the dimensionless electrokinetic width K; that is, the EDL becomes thinner.

Figure 7 depicts the variations of the Nusselt number with zeta potential ζ for different values of the surface heat flux q and the electric field E₀. Both positive and negative values of the wall heat flux are considered, where positive and negative values of q correspond to the case where the fluid is cooled and heated, respectively. In general, increasing the Br decreases the Nusselt number, so the Br decreases with the heat flow flux q. It can be seen that the Brinkman number has a greater influence on the Nusselt number at greater values of heat flux q (see Figure 7a). From Figure 7, we can see that the Nusselt number increases with the positive q because the convective heat transfer coefficient h increases with q. For a given constant wall heat flux q, the Nusselt number is large when E₀ and ζ are small. In contrast to the wall-heating case, in the presence of viscous heating, which occurs at large velocity gradients, q has a greater effect on the Nusselt number. As a result of viscous heating, the bulk temperature is not affected significantly, but the wall temperature is greatly increased, and the Nusselt number will be much higher. Due to the large value of viscous heating, the wall temperature is much higher than the bulk temperature. This behavior of the Nusselt number is accompanied by the occurrence of a singularity in Nusselt number values of the wall heating cases; as a result of changing the sign of dimensionless bulk temperature from negative to positive values (see Figure 7b). Increasing the value of E₀ leads to smaller values of the Nu, indicating that wall temperatures increase rather than bulk temperatures when Joule heating is applied. The reason is that although the distribution of energy generated by Joule heating is uniform throughout the channel cross-section, the energy transferred by convection decreases near the wall, and it equals zero at the wall (see Figure 7c). Moreover, Nu is smaller in the surface charge-dependent slip flow than in the surface charge-independent slip flow. The reason is that the convective heat transfer coefficient h is reduced by the decrease of the velocity on the boundary due to the surface charge effect on the slip.

4. Conclusions

In the present study, the fully developed combined pressure and electroosmotic driven flows in parallel plate micro and nanochannels with surface charge-dependent slip and high zeta potentials have been studied. The classical boundary conditions of surface charge-dependent slip and uniform wall heat flux are taken into account in the analysis. The effects of Joule heating and viscous dissipation are also considered. The electric potential in EDL, velocity, temperature, and Nusselt number have also been obtained analytically. The problem studied here is found to be governed by five parameters: the ratio of the velocity of pressure-driven flow and the velocity of EOF, the electrokinetic width (Debye–Hückel parameter), the zeta potential, Joule heating term, and Brinkman number. The main results of this study can be summarized as follows:

• For a given large zeta potential, the EDL velocity variation is limited to the thin layer that clings to both solid walls.
• Compared with the velocity of the surface charge-independent slip flow, that of the surface charge-dependent slip flow showed a significant reduction.
• The increase in zeta potential causes an increase in the dimensionless temperature in the nanochannel.
• The dimensionless temperature of the surface charge-dependent slip flow is larger than that of the surface charge-independent slip flow.
• For a given zeta potential, the effective Joule heating E_e decreases as the channel height increases; for a given channel height, as the zeta potential ζ increases, the effective Joule heating E_e is increased, especially in the nanochannel.
• Ion density n₀ can inhibit the increase of Joule heat.
• The Nu decreases with Joule heating S and increases with positive heat transfer coefficient q.
• The Nu is large when the applied electric field E₀ and zeta potential ζ are small.
• The Nu in the surface charge-dependent slip flow is smaller than that in the surface charge-independent slip flow.