Study of Heat and Mass Transfer in Electroosmotic Flow of Third Order Fluid through Peristaltic Microchannels

Abstract

An analysis is carried out to evaluate the effects of heat and mass transfer in an electro-osmotic flow of third order fluid via peristaltic pumping. Solutions are derived for small wave number and Peclet number. The emerging non-linear mathematical model is solved analytically and compared numerically by the built-in scheme of working software. The table is inserted for shear stress distribution and a graph for comparison of solution techniques and accuracy of obtained results. The effects of various parameters of interest on pumping, trapping, temperature, heat transfer coefficient, and concentration distribution have been studied graphically. Electro-osmotic exchange of energy and mass has a role in reservoir engineering, chemical industry, and in micro-fabrication technologies.

1. Introduction

In the past few decades, electroosmosis has found rich applications. The first microfluidic devices were developed in the 1980s with the name Micro-Electro-Mechanical Systems (MEMS). These devices have broad application prospects in a variety of fields including biological and medical-related industries, where microfluidic devices are referred to as Lab-on-chip (LOC) devices. LOC can be used in many applications including clinical diagnostics, and biological or chemical contamination. Today, the need for the reliable, fast and well-organized performance of the microfluidic system has created a rigorous demand for small, easy-to-operate and cheap systematic apparatus for drug delivery and DNA analysis.

Microfluidic devices typically use electroosmosis to pump fluid into the entire device. Electroosmosis describes the flow of electrolyte through a channel having a charged boundary due to an applied voltage. Electro-osmosis occurs in many biological, medical and industrial processes such as porous membranes, tubule/canalicular flow, botanical processes, fluid dialysis, transport in human skin and separation techniques. Haung et al. [1] explained the Electroosmotic flow (EOF) in capillaries by monitoring method. Gravesen et al. [2] reviewed the detailed applications of microfluidics. Haswell [3] scrutinized the developmental features of micro-flows on EOF basis. Patankar and Hu [4] predicted the simulation of EOF. Also, Kang et al. [5] examined the EOF in the annulus of capillaries.

The dynamics of peristaltic flow are triggered by its emergence of physiology (intestines, ureters, esophagus, bile ducts, catheters, granules, etc.) and industry (blood pumps, transport of corrosive liquids and hygienic fluid transport). Most of the research on this topic has been conducted on Newtonian fluids. Non-Newtonian fluids cannot be elucidated by a single constitutive relation. Hence, several constitutive relations of these fluids have been proposed. The impact of viscosity in peristalsis for third-grade fluid was explored by Hayat and Abbasi [6]. Noreen [7] proposed the physiological motion of the third order fluid with a magnetic field and mixed convection. Moreover, Noreen [8] considered the peristaltic motion of MHD couple-stress fluid. Prakash et al. [9] worked on peristaltically induced third-grade fluid via the asymmetric conduit. Mallick and Misra [10] discussed Erying–Powell fluid with the electromagnetic field.

Mathematical simulations of peristaltic distribution in microfluidic devices have newly attracted attention. Peristaltic motion can be controlled by adding and resisting external electric fields. One such example is the addition in peristaltic flow by electroosmosis [11]. Tang et al. [12] explored the electroosmosis in the power-law fluid. Hadigol et al. [13] considered the microscopic mixing characteristics of the power-law fluid in electroosmosis. Yavari et al. [14] demonstrated an increase in EOF temperature of bio-fluids (non-Newtonian) via microchannel. Furthermore, Bandopadhyay et al. [15] studied electroosmosis and peristalsis in microchannels. Likewise, Jiang et al. [16] investigated the effect of electroosmosis for Oldroyd-B fluid and wall slip. Jhorar et al. [17] expounded the biomechanical transfer through asymmetric microconduit. Similarly, Tripathi et al. [18] inspected the microvascular blood transport via EOF. Francesko et al. [19] observed the progress of LOC and microfluidics in detail. Javavel et al. [20] illustrated the EOF of pseudoplastic nano liquids through peristaltic pumping. However, the above studies did not examine the electroosmotic effects in peristaltically flowing third order fluid.

Heat transfer in microfluidic devices is of great importance in cancerous tissues destruction, portable kits for diseases diagnostics, the examination of dilution techniques in blood flow and the micro-fabrication technologies. Sadeghi et al. [21] examined the effect of heat exchange in EOF for viscoelastic fluids. Babaie et al. [22] noticed the impact of heat flow for Power-law fluid in the microchannel. Further, Chen et al. [23] analyzed heat exchange for Non-Newtonian fluid suspension in the microchannel. Thermal exchange in blood through capillaries was incorporated by Sinha and Shit [24]. Shit et al. [25] also showed the EOF for heat exchange and MHD. Moreover, in medical operations, it is considered that the wavy walls in electroosmotic flow increase the mass transfer. Therefore, Bhatti et al. [26] described the heat and mass transmission through EDL effects. Also, Reddy et al. [27] revealed the heat transport in peristaltic pumping for Casson fluid through the microchannel. Peristaltic pumping with EDL was studied by Yadav et al. [28]. Narla et al. [29] examined the effect of heat in EOF through time. Moreover, Yang et al. [30] also depicted the transportation of heat in magneto-hydrodynamic EOF.

To our knowledge, there is currently no report on heat and mass transfer in EOF of third-order fluid altered by peristalsis. The present study fills the gap. The main purpose of the present investigation is to analyze the effect of heat transfer and mass transfer for EOF of a third order fluid with viscous dissipation in a microchannel. Note, third order fluid has outstanding significance for shear thickening and shear thinning properties. This investigation has been carried out by using lubrication approximation and Debye Hückel linearization. The non-linear governing equations are solved through perturbation method. Pressure rise and heat transfer coefficient at channel walls have been calculated numerically. Finally, the graphical results for physical quantities are drawn and discussed in detail.

2. Mathematical Model

2.1. Flow Analysis

Consider 2D electroosmotic flow (EOF) of third order fluid in a microfluidic channel having a width $2a$, presented in Figure 1. A peristaltic wave at velocity $c$ propagates along the microchannel wall. The temperature and concentration fields of the lower channel wall are maintained at $T_0$ and $C_0$ whereas upper at $T_1$ and $C_1$, respectively.

Mathematically, microchannel geometry is expressed as:

$$\tilde{Y}=\tilde{h}\left(\tilde{X},\tilde{t}\right) =\left(a+b\sin \left(2\pi \frac{\left(\tilde{X}-c\tilde{t}\right)}{\lambda }\right)\right)$$

(1)

Here $\tilde{X},a, b, c,\lambda$ and $\tilde{t}$ are the axial coordinate, channel’s half width, wave amplitude, wave speed, wavelength and time respectively.

2.2. Governing Equations

In the laboratory frame, governing equations for flow of third order fluid in electrohydrodynamic (EHD) environment are [17]:

Continuity Equation:

$$\frac{\partial \tilde{U}}{\partial \tilde{X}}+\frac{\partial \tilde{V}}{\partial \tilde{Y}}=0,$$

(2)

Equation of motion:

$$\rho \left(\frac{\partial }{\partial \tilde{t}}+\tilde{U}\frac{\partial }{\partial \tilde{X}}+\tilde{V}\frac{\partial }{\partial \tilde{Y}}\right)\tilde{U}=-\frac{\partial \tilde{P}}{\partial \tilde{X}}+\left(\frac{\partial {\tilde{S}}_{\tilde{X}\tilde{X}}}{\partial \tilde{X}}+\frac{\partial {\tilde{S}}_{\tilde{X}\tilde{Y}}}{\partial \tilde{Y}}\right)+{\rho }_e{E}_x,$$

(3)

$$\rho \left(\frac{\partial }{\partial \tilde{t}}+\tilde{U}\frac{\partial }{\partial \tilde{X}}+\tilde{V}\frac{\partial }{\partial \tilde{Y}}\right)\tilde{V}=-\frac{\partial \tilde{P}}{\partial \tilde{Y}}+\left(\frac{\partial {\tilde{S}}_{\tilde{Y}\tilde{X}}}{\partial \tilde{X}}+\frac{\partial {\tilde{S}}_{\tilde{Y}\tilde{Y}}}{\partial \tilde{Y}}\right).$$

(4)

Energy equation with viscous dissipation [9]:

$$\rho c_p\left(\frac{\partial }{\partial \tilde{t}}+\tilde{U}\frac{\partial }{\partial \tilde{X}}+\tilde{V}\frac{\partial }{\partial \tilde{Y}}\right)\tilde{T}=k\left(\frac{{\partial }^2}{\partial {\tilde{X}}^2}+\frac{{\partial }^2}{\partial {\tilde{Y}}^2}\right)\tilde{T}+\left[\frac{\partial \tilde{U}}{\partial \tilde{X}}{\tilde{S}}_{\tilde{X}\tilde{X}}+\frac{\partial \tilde{V}}{\partial \tilde{Y}}{\tilde{S}}_{\tilde{Y}\tilde{Y}}+\left(\frac{\partial \tilde{U}}{\partial \tilde{Y}}+\frac{\partial \tilde{V}}{\partial \tilde{X}}\right){\tilde{S}}_{\tilde{X}\tilde{Y}}\right].$$

(5)

Concentration equation:

$$\left(\frac{\partial }{\partial \tilde{t}}+\tilde{U}\frac{\partial }{\partial \tilde{X}}+\tilde{V}\frac{\partial }{\partial \tilde{Y}}\right)\tilde{C}=D_m\left(\frac{{\partial }^2}{\partial {\tilde{X}}^2}+\frac{{\partial }^2}{\partial {\tilde{Y}}^2}\right)\tilde{C}+\frac{k_t D_m}{T_m}\left(\frac{{\partial }^2}{\partial {\tilde{X}}^2}+\frac{{\partial }^2}{\partial {\tilde{Y}}^2}\right)\tilde{T}.$$

(6)

Here $\tilde{U}$ is the velocity along $\tilde{X}$ direction, $\tilde{V}$ along $\tilde{Y}$ direction, ${\rho }_e$ is electrical charge density, $\tilde{P}$ is pressure, $\rho$ is density of fluid, $E_x$ is axial electric field, $c_p$ is the specific heat, $k$ is thermal conductivity, $k_t$ is ratio of thermal diffusion, $D_m$ is mass diffusivity coefficient, $T_m$ is average temperature, $\tilde{T}$ is temperature and $\tilde{C}$ is concentration field, which are dimensional.

Stress tensor for fluid [6] is defined as

$$\tilde{S}=\left(\mu +{\beta }_3 tr\left({\tilde{A}}_1{}^2\right)\right){\tilde{A}}_1+{\alpha }_1{\tilde{A}}_2+{\alpha }_2{\tilde{A}}_1{}^2+{\beta }_1{\tilde{A}}_3+{\beta }_2\left({\tilde{A}}_1{\tilde{A}}_2+{\tilde{A}}_2{\tilde{A}}_1\right)$$

(7)

where $\mu , {\alpha }_1,{\alpha }_2,{\beta }_1,{\beta }_2$ and ${\beta }_3$ are the material constants. Similarly, Revilin Erickson tensors ${\tilde{A}}_n$ are

$${\tilde{A}}_1=\left(\nabla \tilde{\mathit{V}}\right)+{\left(\nabla \tilde{\mathit{V}}\right)}^{\mathrm{T}}, \mathrm{for} \tilde{\mathit{V}}= [\tilde{U}( \tilde{X}, \tilde{Y},\tilde{t}), \tilde{V}(\tilde{X},\tilde{Y},\tilde{t}), 0]$$

$${\tilde{A}}_n=\frac{\partial {\tilde{A}}_{n-1}}{\partial \tilde{t}}+{\tilde{A}}_{n-1}\left(\nabla \tilde{\mathit{V}}\right)+{\left(\nabla \tilde{\mathit{V}}\right)}^{\mathrm{T}}{\tilde{A}}_{n-1},n>1$$

(8)

The thermodynamic analysis of the model [7] shows that if all the fluid motion is thermodynamically compatible, then these motions must satisfy the Clausius–Duhem inequality. Moreover, it is also supposed that the specific Helmholtz free energy be minimal in equilibrium, then

$$\mu \ge 0, {\alpha }_1\ge 0, \left|{\alpha }_1+{\alpha }_2\right|\le \sqrt{24\mu {\beta }_3} , {\beta }_1={\beta }_2=0, {\beta }_3\ge 0.$$

(9)

In the present investigation, we assume that fluid is thermodynamically compatible, so Equation (7) reduces to the form

$$\tilde{S}=\left(\mu +{\beta }_3 tr\left({\tilde{A}}_1{}^2\right)\right){\tilde{A}}_1+{\alpha }_1{\tilde{A}}_2+{\alpha }_2{\tilde{A}}_1{}^2$$

(10)

2.3. Electrohydrodynamics (EHD)

The Poisson equation [11] in a microchannel is described as:

$${\nabla }^2\tilde{\mathsf{\Phi }}=-\frac{{\rho }_e}{\in }$$

(11)

Here ${\rho }_e,\in$ and $\tilde{\mathsf{\Phi }}$ represents total charge density, dielectric permittivity, and electric potential.

2.4. Potential Distribution

The net charge density ${\rho }_e$ follows the Boltzmann distribution [12], given as

$${\rho }_e=-z_{v}e\left({\tilde{n}}^{-}-{\tilde{n}}^{+}\right),$$

(12)

here the anions $\left({\tilde{n}}^{-}\right)$ and cat-ions $\left({\tilde{n}}^{+}\right)$ are defined through ${\rho }_e$ of the Boltzmann Equation:

$${\tilde{n}}^{\pm }={\tilde{n}}_0{e}^{\left(\pm \frac{e{z}_v}{T_{av}{K}_B}\tilde{\mathsf{\Phi }}\right)},$$

(13)

where $n_0$ represents bulk concentration, $z_v$ the charge balance, $K_B$ the Boltzmann constant, $e$ the electronic charge and $T_{av}$ is the average temperature.

Applying Debye-H$\ddot{u}$ckel linearization approxiamtion [12], Equation (11) becomes

$$\frac{d^2\mathsf{\Phi }}{d{y}^2}=m_e{}^2\mathsf{\Phi }.$$

(14)

Here $m_e$ is the electroosmotic parameter. The analytical solution of above Equation (14) subject to boundary conditions $\frac{\partial \mathsf{\Phi }}{\partial y}=0,$ at $y=0$ and $\mathsf{\Phi }=1,$ at $y=h\left(x\right)$ is obtained as;

$$\mathsf{\Phi }\left(y\right)=\frac{\cosh \left(m_{e}y\right)}{\cosh \left(m_{e}h\right)}.$$

(15)

2.5. Non-Dimensionlization, Lubrication Approach and Boundary Conditions

In the stationary frame $(\tilde{X},\tilde{Y})$, the motion is unsteady due to the moving boundary. But, if viewed in a moving frame $\left(\tilde{x},\tilde{y}\right)$, it can be considered as steady due to static boundary. Therefore, the translational transformation between two frames is [7]:

$$\begin{gathered} \tilde{x}=\tilde{X}-c\widetilde{t,} \tilde{y}=\tilde{Y}, \tilde{u}(\tilde{x},\tilde{y})=\tilde{U}(\tilde{X},\tilde{Y},\tilde{t})-c, \tilde{v}(\tilde{x},\tilde{y})=\tilde{V}(\tilde{X},\tilde{Y},\tilde{t}), \\ \tilde{p}(\tilde{x},\tilde{y})=\tilde{P}(\tilde{X},\tilde{Y},\tilde{t}), T(\tilde{x},\tilde{y})=\tilde{T}(\tilde{X},\tilde{Y},\tilde{t}) \end{gathered}$$

(16)

The above transformations is used in Equations (2)–(6) and then introducing the non-dimensional variables:

$$\begin{gathered} x=\frac{2\pi \tilde{x}}{\lambda }, y=\frac{\tilde{y}}{a}, t=\frac{2\pi c\tilde{t}}{\lambda }, u=\frac{\tilde{u}}{c}, v=\frac{\tilde{v}}{c}, p=\frac{2\pi a^2}{\lambda c\mu }\tilde{p},\alpha =\frac{2\pi a}{\lambda }, h=\frac{\tilde{h}}{a},\mathsf{\Phi }=\frac{\tilde{\mathsf{\Phi }}}{a}, \\ {P}_e=\frac{c\lambda }{D},m_e=\frac{a}{{\lambda }_D},n=\frac{\tilde{n}}{n_0},\mathsf{\epsilon }=\frac{b}{a}, \mathsf{\theta }=\frac{\tilde{T}-T_0}{T_0},\Omega =\frac{\tilde{C}-C_0}{C_0},{\lambda }_D=\frac{1}{e{z}_v}\sqrt{\frac{T_0\in K_B}{2n_0}}, \\ E_c=\frac{c^2}{c_p{T}_0},\mathsf{\Phi }=\frac{e{z}_v}{T_{av}{K}_B}\tilde{\mathsf{\Phi }}, U_{hs}=-\frac{E_x\in }{c\mu }, S=\frac{a}{c\mu }\tilde{S}, {\lambda }_1=\frac{{\alpha }_{1}c}{a\mu },{\lambda }_2=\frac{{\alpha }_{2}c}{a\mu },{\gamma }_1=\frac{{\beta }_1{c}^2}{a^2\mu }, \\ {S}_c=\frac{\nu }{D_m}, {\gamma }_2=\frac{{\beta }_2{c}^2}{a^2\mu }, {\gamma }_3=\frac{{\beta }_3{c}^2}{a^2\mu }, \mathsf{\psi }=\frac{\tilde{\mathsf{\psi }}}{ca}, u=\frac{\partial \mathsf{\psi }}{\partial y},v=-\alpha \frac{\partial \mathsf{\psi }}{\partial x},R_e=\frac{\rho ca }{\mu }, \\ {P}_r=\frac{\mu c_p}{k}, B_r=E_c.P_r, \Gamma ={\gamma }_3+{\gamma }_2, \mathrm{f}=\frac{\tilde{\mathrm{q}}}{ca},\mathsf{\Theta }=\frac{\tilde{\mathrm{Q}}}{ca},S_r=\frac{T_0 k_t D_m}{\nu T_m C_0}{.}_{ } \end{gathered}$$

(17)

Here, $x$ and $y$ are the axial and transverse coordinate, $m_e$ is electroosmotic parameter, $p$ is dimensionless pressure, ${\lambda }_D$ is Debye length, $\alpha$ is peristaltic wave number, $\mathsf{\Theta }$ is dimensionless volume flow rate, $\mathsf{\theta }$ is dimensionless temperature, $\Omega$ is dimensionless concentration, $U_{hs}$ is Helmholtz–Smoluchowski velocity, $\mathsf{\psi }$ is non-dimensional stream function, $P_e$ is ionic Peclet number, $B_r$ is Brinkman, $R_e$ is Reynolds number, $P_r$ is Prandtl, $\Gamma$ is Deborah number, $\mathsf{\epsilon }$ is amplitude ratio, $S$ is non-dimensional shear stress, $S_r$ is non-dimensional Soret number, $S_c$ is Schmidt number and $E_c$ is Eckert number.

The relation of velocity and stream function $\mathsf{\psi }$ is defined as

$$u=\frac{\partial \mathsf{\psi }}{\partial y}, v=-\alpha \frac{\partial \mathsf{\psi }}{\partial x}.$$

(18)

After introducing Equations (17) and (18), Equation (2) is satisfied identically and Equations (3)–(6) becomes

$$R_e\alpha \left(\frac{\partial \mathsf{\psi }}{\partial y}\frac{\partial }{\partial x}-\frac{\partial \mathsf{\psi }}{\partial x}\frac{\partial }{\partial y}\right)\frac{\partial \mathsf{\psi }}{\partial y}=-\frac{\partial p}{\partial x}+\left(\alpha \frac{\partial S_{xx}}{\partial x}+\frac{\partial S_{xy}}{\partial y}\right)+m_e{}^2{U}_{hs}\mathsf{\Phi },$$

(19)

$$-R_e{\alpha }^3\left(\frac{\partial \mathsf{\psi }}{\partial y}\frac{\partial }{\partial x}-\frac{\partial \mathsf{\psi }}{\partial x}\frac{\partial }{\partial y}\right)\frac{\partial \mathsf{\psi }}{\partial x}=-\frac{\partial p}{\partial y}+\alpha \left(\alpha \frac{\partial S_{xy}}{\partial x}+\frac{\partial S_{yy}}{\partial y}\right),$$

(20)

$$R_e\alpha \left(\frac{\partial \mathsf{\psi }}{\partial y}\frac{\partial }{\partial x}-\frac{\partial \mathsf{\psi }}{\partial x}\frac{\partial }{\partial y}\right)\mathsf{\theta }=\frac{1}{P_r}\left({\alpha }^2\frac{{\partial }^2}{\partial x^2}+\frac{{\partial }^2}{\partial y^2}\right)\mathsf{\theta }+E_c\left[{\alpha }^2\frac{\partial }{\partial x}\left(\frac{\partial \mathsf{\psi }}{\partial y}\right)S_{xx}-\alpha \frac{\partial }{\partial y}\left(\frac{\partial \mathsf{\psi }}{\partial x}\right)S_{yy}+\left(\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}-{\alpha }^2\frac{{\partial }^2\mathsf{\psi }}{\partial x^2}\right)S_{xy}\right],$$

(21)

$$R_e\alpha \left(\frac{\partial \mathsf{\psi }}{\partial y}\frac{\partial }{\partial x}-\frac{\partial \mathsf{\psi }}{\partial x}\frac{\partial }{\partial y}\right)\Omega =\frac{1}{S_c}\left({\alpha }^2\frac{{\partial }^2}{\partial x^2}+\frac{{\partial }^2}{\partial y^2}\right)\Omega +S_r\left({\alpha }^2\frac{{\partial }^2}{\partial x^2}+\frac{{\partial }^2}{\partial y^2}\right)\mathsf{\theta }$$

(22)

Now subject to long-wavelength assumption [10] on Equations (19)–(22), ignoring higher powers of $\alpha$, we get;

$$\frac{\partial p}{\partial x}=\frac{\partial }{\partial y}\left[\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}+2\Gamma {\left(\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}\right)}^3\right]+m_e{}^2{U}_{hs}\mathsf{\Phi },$$

(23)

$$\frac{\partial p}{\partial y}=0$$

(24)

$$\frac{{\partial }^2\mathsf{\theta }}{\partial y^2}+B_r\left[{\left(\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}\right)}^2+2\Gamma {\left(\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}\right)}^4\right]=0$$

(25)

$$\frac{{\partial }^2\Omega }{\partial y^2}+S_c.S_r\left(\frac{{\partial }^2\mathsf{\theta }}{\partial y^2}\right)=0$$

(26)

where $S_{xy}=\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}+2\Gamma {\left(\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}\right)}^3$.

Eliminating pressure from Equations (23) and (24), we get

$$\frac{{\partial }^2}{\partial y^2}\left[\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}+2\Gamma {\left(\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}\right)}^3\right]+m_e{}^2{U}_{HS}\frac{\partial \mathsf{\Phi }}{\partial y}=0$$

(27)

The dimensionless boundary conditions are:

$$\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}=0, \mathsf{\psi }=0, \frac{\partial \mathsf{\theta }}{\partial y}=0, \frac{\partial \Omega }{\partial y}=0, \mathrm{at} y=0{,}_{ }$$

(28)

$$\frac{\partial \mathsf{\psi }}{\partial y}=-1, \mathsf{\psi }=\mathrm{f}, \mathsf{\theta }=0, \Omega =0, \mathrm{at} y=h\left(x\right)=1+\mathsf{\epsilon }\sin x_{ }$$

(29)

here $\mathrm{f}$ is dimensionless mean flow rate in moving frame and $\mathsf{\Theta }$ is dimensionless mean flow rate in fixed frame. Furthermore, $\mathrm{f}$ in the moving frame is connected to the $\mathsf{\Theta }$ in the fixed frame as:

$$\mathsf{\Theta }=\mathrm{f}+1.$$

(30)

here $\mathrm{f}={{\displaystyle \int }}_0^{h\left(x\right)}\frac{\partial \mathsf{\psi }}{\partial y}dy=\mathsf{\psi }\left(h\left(x\right)\right)-\mathsf{\psi }\left(0\right)$.

3. Solution Methodology

3.1. Series Solutions

The governing differential system given in Equations (23)–(27) subject to boundary conditions (28,29) consists of highly non-linear and coupled equations. Exact solution of the governing system of equations is impossible. Therefore, analytical solution of the system is presented through the perturbation technique by taking Deborah number $\Gamma$ as a perturbation parameter. So, we expand $\mathsf{\psi }$, $\mathsf{\theta }$, $\Omega ,$. $dp/dx$ and $\mathrm{f}$ as in [8]:

$$\mathsf{\psi }={\mathsf{\psi }}_0+\Gamma {\mathsf{\psi }}_1+O{\left(\Gamma \right)}^2$$

(31)

$$\mathsf{\theta }={\mathsf{\theta }}_0+\Gamma {\mathsf{\theta }}_1+O{\left(\Gamma \right)}^2$$

(32)

$$\Omega ={\Omega }_0+\Gamma {\Omega }_1+O{\left(\Gamma \right)}^2$$

(33)

$$dp/dx=d{p}_0/dx+\Gamma d{p}_1/dx+O{\left(\Gamma \right)}^2$$

(34)

$$\mathrm{f}={ \mathrm{f}}_0+\Gamma {\mathrm{f}}_1+{\left(\Gamma \right)}^2$$

(35)

3.2. Zero Order System ${\left(\Gamma \right)}^0$

$$\frac{{\partial }^4{\mathsf{\psi }}_0}{\partial y^4}+m_e{}^2{U}_{HS}\frac{\partial \mathsf{\Phi }}{\partial y}=0,$$

(36)

$$\frac{{\partial }^2{\mathsf{\theta }}_0}{\partial y^2}+B_r{\left(\frac{{\partial }^2{\mathsf{\psi }}_0}{\partial y^2}\right)}^2=0$$

(37)

$$\frac{{\partial }^2{\Omega }_0}{\partial y^2}+S_c.S_r\left(\frac{{\partial }^2{\mathsf{\theta }}_0}{\partial y^2}\right)=0$$

(38)

$$\frac{\partial p_0}{\partial x}=\frac{{\partial }^3{\mathsf{\Psi }}_0}{\partial y^3}+m_e{}^2{U}_{HS}\mathsf{\Phi },$$

(39)

subject to

$$\begin{gathered} {\mathsf{\psi }}_0=0, \frac{{\partial }^2{\mathsf{\psi }}_0}{\partial y^2}=0, \frac{\partial {\mathsf{\theta }}_0}{\partial y}=0, \frac{\partial {\Omega }_0}{\partial y}=0, \mathrm{at} y=0 \\ {\mathsf{\psi }}_0={ \mathrm{f}}_0, \frac{\partial {\mathsf{\psi }}_0}{\partial y}=-1, {\mathsf{\theta }}_0=0, {\Omega }_0=0, \mathrm{at} y=h \end{gathered}$$

(40)

3.3. First Order System ${\left(\Gamma \right)}^1$

$$\frac{{\partial }^4{\mathsf{\psi }}_1}{\partial y^4}+6{\left(\frac{{\partial }^2{\mathsf{\psi }}_0}{\partial y^2}\right)}^2\frac{{\partial }^4{\mathsf{\psi }}_0}{\partial y^4}+12{\left(\frac{{\partial }^3{\mathsf{\psi }}_0}{\partial y^3}\right)}^2\frac{{\partial }^2{\mathsf{\psi }}_0}{\partial y^2}=0,$$

(41)

$$\frac{{\partial }^2{\mathsf{\theta }}_1}{\partial y^2}+2B_r\left[\left(\frac{{\partial }^2{\mathsf{\psi }}_0}{\partial y^2}\right)\left(\frac{{\partial }^2{\mathsf{\psi }}_1}{\partial y^2}\right)+{\left(\frac{{\partial }^2{\mathsf{\psi }}_0}{\partial y^2}\right)}^4\right]=0$$

(42)

$$\frac{{\partial }^2{\Omega }_1}{\partial y^2}+S_c.S_r\left(\frac{{\partial }^2{\mathsf{\theta }}_1}{\partial y^2}\right)=0$$

(43)

$$\frac{\partial p_1}{\partial x}=\frac{{\partial }^3{\mathsf{\psi }}_1}{\partial y^3}+6{\left(\frac{{\partial }^2{\mathsf{\psi }}_0}{\partial y^2}\right)}^2\frac{{\partial }^3{\mathsf{\psi }}_0}{\partial y^3},$$

(44)

subject to

$$\begin{gathered} {\mathsf{\psi }}_1=0, \frac{{\partial }^2{\mathsf{\psi }}_1}{\partial y^2}=0, \frac{\partial {\mathsf{\theta }}_1}{\partial y}=0, \frac{\partial {\Omega }_1}{\partial y}=0, \mathrm{at} y=0, \\ {\mathsf{\psi }}_1={ \mathrm{f}}_1, \frac{\partial {\mathsf{\psi }}_1}{\partial y}=0, {\mathsf{\theta }}_1=0, {\Omega }_1=0, \mathrm{at} y=h \end{gathered}$$

(45)

The solution for zero and first order system is defined in the appendix. By substituting the zero and first order solution in the Equations (32)–(35) we get the final first order solution for $\mathsf{\psi }$, $\mathsf{\theta }$, $\Omega$ and $dp/dx$.

The dimensionless pressure rises per wavelength $\Delta P_{\lambda }$ is as follows:

$$\Delta P_{\lambda }={{\displaystyle \int }}_0^1\frac{\partial p}{\partial x}dx.$$

(46)

The heat transfer coefficient $Z$ is defined as:

$$Z=h_x\frac{\partial \mathsf{\theta }}{\partial y}{|}_{y=h\left(x\right)}.$$

(47)

The shear stress $S_{xy}$ at the wall $y=h\left(x\right)$ is defined as

$$S_{xy}=\left(\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}+2\Gamma {\left(\frac{{\partial }^2\mathsf{\psi }}{\partial y^2}\right)}^3\right){|}_{y=h\left(x\right)}.$$

(48)

4. Computational Discussion of Results

This section deals with the impacts of different parameters (i.e., Deborah number $\Gamma$, electroosmotic parameter $m_e$, Helmholtz–Smoluchowski/maximum electroosmotic velocity $U_{hs}$, Brinkman number $B_r$, volume flow rate $\mathsf{\Theta }$, Soret number $S_r$ and Schmidt number $S_c$) on the pressure rise $\Delta P_{\lambda }$, pressure gradient $dp/dx$, velocity distribution $u$, temperature distribution $\mathsf{\theta }$, concentration distribution $\Omega$, heat transfer coefficient $Z$ and streamlines $\mathsf{\psi }$. Figure 2, Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9, Figure 10 and Figure 11 are plotted to serve the purpose.

4.1. Comparative Analysis

Figure 2 reveals the overlapping comparison between the analytical and numerical solution of velocity. Note that the analytical solution is obtained through perturbation technique and numerical solution is from built-in scheme NDsolve of working software. It is concluded that our results through perturbation (i.e., truncated up to first order of third-order fluid parameter) are in good agreement with the numerical scheme.

4.2. Flow Characteristics

The axial velocity deserves to provide the main features of flow behavior in the microchannel for microfluidic applications. Figure 3a–c is portrayed to inspect the evolutions in the velocity profile across the microchannel for various values of Deborah number $\Gamma$, electroosmotic parameter $m_e$ and maximum electroosmotic velocity $U_{hs}$. Figure 3a shows the variation of axial velocity for various values of $\Gamma$. It is demonstrated that the magnitude of velocity $u$ decreases as $\Gamma$ increases at $y=0$. The reason behind this trend is EDL (electric double layer). It means flow of fluid resists in the presence of EDL. But opposite behavior is witnessed near the channel wall $y=h\left(x\right)$. Figure 3b predicts that $m_e$ boosts the magnitude of $u$ i.e., magnitude of $u$ increases as $m_e$ increases at the central region of the channel and reduces near the walls of the channel. Since $m_e$ is the ratio of the channel height to the Debye length ${\lambda }_D$, it specifies that ${\lambda }_D$ is inversely proportional to EDL. Hence, more fluid flows in the central region $y=0$. Figure 3c depicts that $u$ decent as $U_{hs}$ ascents in $y=0$, whereas opposite movement is noticed at $y=h\left(x\right)$. As $U_{hs}$ physically concludes that velocity of fluid decreases when thickness of EDL increases. Therefore, fluid flow decreases in EDL presence.

4.3. Pumping Characteristics

Peristalsis is characterized by pushing the fluid forward against the pressure rise $\Delta P_{\lambda }$. The characteristic of pumping can be analyzed in the form of $\Delta P_{\lambda }$ versus $\mathsf{\Theta },$ shown in Figure 4a–c Pumping action divides the whole region into the four segments: pumping region (adverse pressure gradient $\Delta P_{\lambda }>0$, Positive pumping $\mathsf{\Theta }>0$, backward/retrograde pumping ($\Delta P_{\lambda }>0$,$\mathsf{\Theta }<0$), augmented pumping (favorable pressure gradient $\Delta P_{\lambda }<0$, Positive pumping $\mathsf{\Theta }>0$) and free pumping ($\Delta P_{\lambda }=0$). Figure 4a portrays the $\Delta P_{\lambda }$ against $\mathsf{\Theta }$ for different values of $\Gamma$. It is perceived that pumping increases as $\Gamma$ increases in pumping region ($\Delta P_{\lambda }>0$). But in co-pumping ($\Delta P_{\lambda }<0$), pumping decreases as $\Gamma$ increases. It means that pumping rate is high for third order fluid under parameter $\Gamma$ in comparison to viscous fluid $\Gamma =0$. For free pumping ($\Delta P_{\lambda }=0$), curves coincide i.e., there is no difference between third-grade fluid and Newtonian fluid within the domain of −0.4 to 0.4. Figure 4b describes that $\Delta P_{\lambda }$ decreases as $m_e$ increases. This means additional pressure is needed to drive the third-grade fluid for maximum thickness of EDL. It physically declares that formation of EDL on charged surface opposes the third-grade fluid flow. Furthermore, Figure 4c illustrates $\Delta P_{\lambda }$ decreases as $U_{hs}$ increases. Here pumping is controlled through an electric field. Note that there is a linear relationship between $\Delta P_{\lambda }$ against $\mathsf{\Theta },$ because of $U_{hs}$.

The distribution of pressure gradient $dp/dx$ is plotted through Figure 5a–d for Deborah number $\Gamma$, electroosmotic parameter $m_e$, maximum electroosmotic velocity $U_{hs}$ and volume flow rate $\mathsf{\Theta }$. Figure 5a portrays that by increasing the $\Gamma ,$ magnitude of pressure gradient $dp/dx$ increase. $\Gamma$ being a physical parameter showing the Non-Newtonian nature, one can observe pressure gradients are higher for third order in comparison to Newtonian fluids. Figure 5b illustrates that by increasing the values of $\mathsf{\Theta },$. $dp/dx$ decreases. Because there is an inverse relationship between $dp/dx$ and $\mathsf{\Theta }$. Figure 5c highlights that by raising the values of $m_e,$. $dp/dx$ decreases. It means that EDL presence in charged surface resists the flow of third grade fluid. Similarly, it is important that the characteristics of pumping can be modified by EDL mechanism and process of pumping can be systematized by thinning and thickening EDL width. Also, Figure 5d presents that the magnitude of $dp/dx$ increases as $U_{hs}$ increases.

4.4. Trapping Characteristics

Trapping for different values of Deborah number $\Gamma$, electroosmotic parameter $m_e$ and maximum electroosmotic velocity $U_{hs}$ are shown in Figure 6, Figure 7 and Figure 8. Figure 6a–d reveal a streamline structure for various values of $\Gamma$. The size of the enclosed bolus decreases as $\Gamma$ increases. It means that size of trapped bolus is strongly affected by changing the fluid. Similarly, Figure 7a–d visualize that the number of enclosed bolus increases with rise in electroosmotic parameter $m_e$. It scrutinizes that as the $m_e$ increases, the enclosed bolus strongly appears in EDL. It means that more fluid can be trapped in the presence of more electric field. Figure 8a–d illustrate streamline makeup for $U_{hs}$. It demonstrates that the accumulation of streamline reduces with rise in $U_{hs}$. Thus, stronger $U_{hs}$ means stronger the external electric field, the number of enclosed significantly decreases for large values of bolus $U_{hs}$.

4.5. Temperature Characteristics

In the present subsection, we examined the effect of Deborah number $\Gamma$, electroosmotic parameter $m_e$, maximum electroosmotic velocity $U_{hs}$ and Brinkman number $B_r$ on heat transfer characteristics. Figure 9a reveals the various values of $\Gamma$, for temperature distribution $\mathsf{\theta }.$ It examines that $\mathsf{\theta }$ enhances for more values of $\Gamma .$ Also, Figure 9b displays that $\mathsf{\theta }$ increases by increasing $m_e$. Therefore, the decrease in EDL produces a rise in $\mathsf{\theta }$ remarkably. Figure 9c, disclose the effects of maximum electroosmotic velocity $U_{hs}$ on $\mathsf{\theta }$. It signifies that for more values of $U_{hs}$, $\mathsf{\theta }$ declines. Figure 9d exemplifies the variation of $B_r$ on $\mathsf{\theta }$ i.e., $\mathsf{\theta }$ increases for more values of $B_r$. The main fact in $B_r$ is the viscosity, which creates resistance. This resistance is responsible for the collision of fluid particles and the resulting collision causes an increase in temperature.

4.6. Concentration Characteristics

Here, we observed the impact of Deborah number $\Gamma$, electroosmotic parameter $m_e$, maximum electroosmotic velocity $U_{hs}$, Brinkman number $B_r$, Schmidt number $S_c$ and Soret number $S_r$ on mass transfer characteristics. Figure 10a reveals the various values of $\Gamma$, for temperature distribution $\mathsf{\theta }.$ It examines that $\Omega$ decreases for more values of $\Gamma .$ Figure 10b displays that $\Omega$ decreases by increasing $m_e$. Figure 10c discloses the effects of $U_{hs}$ on $\Omega$. It signifies that for more values of $U_{hs}$, $\mathsf{\theta }$ rises. Figure 10d exemplifies the variation of $B_r$ on $\Omega$ i.e., $\Omega$ decreases for more values of $B_r$. Figure 10e discloses the effects of $S_c$ on $\Omega$. It implies that $\Omega$ decreases for more values of $S_c.$ Also, Figure 10f shows that $\Omega$ decreases for more values of $S_r.$

4.7. Heat Transfer Coefficient

Figure 11a–d demonstrates the changes in heat transfer coefficient $Z$ for various values of Deborah number $\Gamma$, electroosmotic parameter $m_e$, maximum electroosmotic velocity $U_{hs}$, and Brinkman number $B_r$. This reveals that the various values of $\Gamma$, $B_r$ and $m_e$, the heat transfer coefficient increases. On contrary $Z$ is decreased for different values of $U_{hs}$.

4.8. Shear Stress Distribution

Distribution of shear stress provides very useful information about the nature of dissipation at walls.

Table 1. Numerical values of shear stress at the wall $y=1+\mathsf{\epsilon }\sin x$ for different values of Deborah number $\Gamma$, Helmholtz–Smoluchowski velocity $U_{hs}$ and electroosmotic parameter $m_e$.

(a)$\mathsf{\epsilon }\mathbf{=}\mathbf{0.385}\mathbf{,}\mathbf{\Theta }\mathbf{=}\mathbf{1.0}\mathbf{,}{\mathit{m}}_{\mathit{e}}\mathbf{=}\mathbf{2}\mathbf{,}{\mathit{U}}_{\mathit{h}\mathit{s}}\mathbf{=}\mathbf{-}\mathbf{1.0}\mathbf{.}$
$\mathit{x}$	$\mathit{\Gamma }$
	0.0	0.2	0.4
0.1	−2.4917	−3.0742	−0.8223
0.2	−2.3669	−3.8790	−2.2711
0.3	−2.6497	−4.5689	−4.3309
(b)$\mathsf{\epsilon }\mathbf{=}\mathbf{0.385}\mathbf{,}\mathbf{\Theta }\mathbf{=}\mathbf{1.0}\mathbf{,}{\mathit{m}}_{\mathit{e}}\mathbf{=}\mathbf{2}\mathbf{,}\mathit{\Gamma }\mathbf{=}\mathbf{0.1}\mathbf{.}$
$\mathit{x}$	${\mathit{U}}_{\mathit{h}\mathit{s}}$
	−1.0	0.0	1.0
0.1	−3.6869	−1.1370	0.5727
0.2	−3.7974	−1.2903	0.4311
0.3	−3.8395	−1.4188	0.3188
(c)$\mathsf{\epsilon }\mathbf{=}\mathbf{0.385}\mathbf{,}\Theta \mathbf{=}\mathbf{1.0}\mathbf{,}\mathit{\Gamma }\mathbf{=}\mathbf{0.1}\mathbf{,}{\mathit{U}}_{\mathit{h}\mathit{s}}\mathbf{=}\mathbf{-}\mathbf{1.0}$
$\mathit{x}$	${\mathit{m}}_{\mathit{e}}$
	2	4	6
0.1	−3.6869	−10.778	−6.3530
0.2	−3.7974	−8.7149	−4.7301
0.3	−3.8395	−7.0505	−3.6590

a–c show the distribution of axial shear stress for different values of Deborah number $\Gamma$, maximum electroosmotic velocity $U_{hs}$ and electroosmotic parameter $m_e$. Table 1a explains that the magnitude of shear stress decreases in the axial direction while it decreases (near the channel wall) then increases (away from channel wall) for different values of $\Gamma .$ Here, it is obvious that the third-grade fluid act as a shear thinning and shear thickening fluid. Similarly, Table 1b displays the variation of shear stress for $U_{hs}$. It is found that the axial distribution of shear stress enhances for increasing values of $U_{hs}$. Table 1c presents the shear stress distribution for different values of $m_e.$ It is concluded that the same behavior is observed as of $\Gamma$. Therefore, we say that the presence of electroosmotic parameter strongly infers the characteristics of third order fluid parameter.

5. Concluding Remarks

From current analysis, we conclude that:

• The axial velocity of third-grade fluid in the microchannel is enhanced due to an increase in electroosmotic parameter.
• For higher values of Schmidt or Soret number, concentration profile decreases.
• The magnitude of pressure rise decreases in the pumping region with the increase of electroosmotic parameter.
• Pressure gradient is more for EOF of third order fluid in comparison to peristaltically flowing fluid.
• Axial distribution of shear stress enhances for increasing values of maximum electroosmotic velocity.
• Volume of trapped bolus is strongly dependent upon the electroosmotic parameter whereas it is suppressed for Deborah number and maximum electroosmotic velocity.
• Temperature distribution in microchannel is crucial for electroosmotic parameter, Deborah number and Brinkman number.

Electro-osmotic exchange of energy and mass has a role in reservoir engineering. The outcomes of the present model may be used by engineers in chemical industry and in micro-fabrication technologies. In the future, this model accompanying porosity will help in understanding the hydrodynamics of rheological fluids in a typical channel employed in a LOC system. This course requires promotion, being of significance in factual life applications.