Effectiveness of flow obstructions in enhancing electro-osmotic flow

In this paper the influence of obstructions on microchannel electro-osmotic flow is investigated for the first time. To carry out such a study, regular obstructions are introduced into microchannels and flow rates are numerically calculated. The effect of channel width on flow rates is analysed on both free and obstructed channels. The solid material considered for channel walls and obstructions is silicon, and the electrolyte is deionised water. The parameters studied include channel width, obstruction size and effective porosity of the channel. The effective porosity is varied between 0.4 and 0.8 depending on other chosen parameters. The results clearly demonstrate that, under the analysed conditions, introduction of obstructions into channels wider than 100μm\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$100\,\upmu \hbox {m}$$\end{document} enhances the flow rate induced by electro-osmosis.

, while others have used a generalised model for porous media flow and added a source term in the momentum equation, depending on the charge density of porous medium and the applied electrical field (Scales 2004;Tang et al. 2010). Although it has been found that the main driving force is due to the charged particles rather than the channel walls (Wang et al. 2006), it appears that the internal potential equation has not been appropriately modified to take into account the charge of solid particles, except for boundary conditions (Tang et al. 2010). To consider the charge of both solid particles and channel walls, Kang et al. (2005) split the velocity into two components and then coupled them to obtain the overall macroscopic EOF velocity. The first component was derived as per the fluid flow in standard channels, by assimilating the porous medium to an assembly of parallel tortuous cylinders. The second component was obtained by applying the Brinkman extension of the Darcy equation, in which the inertia terms were neglected because of the low Reynolds number. The dimensionless Darcy velocity was found to increase with the particles size, the applied electric field and the difference between zeta potential of particles and channel walls, and it was also found to decrease with increase in channel width.
Many authors have focused on EO porous pumps and found that the thermodynamic efficiency significantly increases with the addition of a porous medium in a channel, as much higher pumping pressures can be generated (Wang et al. 2006).
In general, the velocity has been found to increase with the increase in diameter of solid particles or pores (Berrouche et al. 2009;Chai et al. 2007;Chen et al. 2008;Kang et al. 2005;Tang et al. 2010;Wang and Chen 2007;Yao et al. 2006) and porosity (Chai et al. 2007;Tang et al. 2010;Wang and Chen 2007). Although EOF can be enhanced by increasing the contact surface between the electrolyte and the boundaries of charged particles, excessive reduction in void space between the particles (Kang et al. 2005) can substantially increase the viscous losses (Berrouche et al. 2009). Such losses can easily compromise the pressure generated by surface charge (Chai et al. 2007). The effect of the distance between solid particles on EOF has also been analysed in the past. By considering only the charge of solid particles, the maximum velocity was achieved for a distance between particles equal to five times that of the Debye length (see Sect. 2.2) (Li et al. 2013a). It appears that the EO driving force decreases when the distance between the particles is different from this value (Li et al. 2013a). When the gap between the particles is too small, the electric double layer of the particles overlaps and as a consequence the velocity decreases several orders of magnitude. An attempt to optimise the porosity has found that as the pore number density increases and the pore diameter decreases, both flow rate and pressure performance of porous membrane pumps increase (Yao et al. 2006). Cheema et al. (2013) found that EO pumping can be improved when the layers adjacent to the solid wall of the pump present higher porosity with respect to the central region. In this case, the velocity near the walls is larger and therefore variation in porosity near the wall cannot be neglected (Chai et al. 2007). Chen et al. (2014) investigated the influence of the structure of porous media on EO permeability, that relates the average fluid velocity to the external electric field, and found that EO permeability grows monotonically with increasing porosity.
As mentioned before, EOF in porous media has been modelled by using a microscopic approach in many studies. It was found that the geometry used to approximate the morphology of porous media significantly affects the prediction of the volumetric flow rate that was observed to be greater in cylindrical and annular capillaries than in rectangular channels (Pascal et al. 2012). Morphology of porous media was also investigated by Wang (2012), who analysed the difference in terms of EOF permeability obtained by using granular, fibrous or network structures. He showed that at low porosity the network configuration enhanced EOF permeability, due to its highest surface/volume ratio, whereas at higher values of porosity the granular structure performed better, for its lower resistance to flow.
From the analysis of the available scientific literature on the subject of EOF in porous media, it appears that the work on that subject is quite heterogeneous in nature (Di Fraia et al. 2017). Although some important attempts have been made in the recent past, many of the fundamental aspects are still unclear. The two important issues that need a better understanding are: method of modelling EOF past obstructions and flow enhancement as a result of introducing obstructions confined in a microchannel. By focusing on these aspects the effectiveness of using flow obstructions in microchannels to enhance EOF should be determined with certainty.
The aim of this paper is to understand under which conditions the use of obstructions in microchannels can be useful to increase EO flow rate. For this reason, in the present work a well-defined porous medium represented by regular flow obstructions is studied by using a microscopic approach that provides details of flow structure at particle level (Massarotti et al. 2003). The interaction between obstruction boundaries and side walls, which affects the internal potential distribution, and consequently EOF, has been analysed. Fluid flow is modelled by using Navier-Stokes equations for incompressible flow that are solved by using a fully explicit artificial compressibility-based CBS (characteristic-based split) scheme (Nithiarasu 2003;Nithiarasu et al. 2016) through the channel. The momentum equation has been modified to include the electro-kinetic effects responsible for EOF into a source term. The effect of obstructions on the internal potential distribution is taken into account by introducing a novel reference length. To analyse the effectiveness of introducing porous media to enhance EOF, simulations have been carried out on microchannels with and without obstructions and the results are compared in terms of internal potential distribution, fluid velocity and flow rate. Furthermore, the model has been used to carry out a sensitivity analysis of these quantities on the channel width and on the topology of the obstructions.
The paper is organised into the following sections. In Sect. 2 the governing equations are briefly presented (for the sake of completeness the solution algorithm is described in "Appendix"). In Sect. 3 the numerical results are reported: firstly the results for fluid flow without obstructions are compared to their available analytical solutions; then, the model is used to assess the effect of channel width, particle size and porosity of the channel with obstructions on EOF. Finally, in Sect. 4, some conclusions are drawn.

Mathematical model and solution procedure
The electric field that generates electro-osmotic flow (EOF) is induced by the interaction between an external applied potential and the electric double layer (EDL), as shown in Fig. 1.
The electro-kinetic forces responsible for EOF are modelled through Laplace equation governing the external applied potential and Poisson-Boltzmann equation governing the EDL potential. These equations defining the electrical field are decoupled and solved separately from the Navier-Stokes equations, and their effect on the flow is taken into account through a source term of the momentum equation (Patankar and Hu 1998;Yang and Li 1998). Fluid flow is analysed through the Navier-Stokes equations that are temporally discretised by using the characteristic-based split (CBS) algorithm (Nithiarasu 2003;Nithiarasu et al. 2016), while the Galerkin approximation is used for spatial discretisation.

Governing equations
The external potential φ is governed by the Laplace equation of the type: where σ is the fluid electrical conductivity. The external electric field, E x , and the external electric potential, φ, are related via The electric double layer (EDL) potential, ψ, is described by the Poisson-Boltzmann equation, as follows: where ε is the dielectric constant of the electrolyte, ε 0 is the permittivity of vacuum, and ρ E is the net charge density.
The equilibrium Boltzmann distribution equation can be used to predict the ionic number concentration in the case of fully developed flow with small Peclet numbers (Yang et al. 2001) and in the absence of overlapping of EDLs (Qu and Li 2000). The ions in the solution are assumed to be equal but opposite in charge, and they are related to their energy (zeψ) as where n + and n − represent the number of positive and negative ions, respectively, n 0 is the ionic number concentration in the bulk solution, z is the valance of the ions, e is the elementary charge, k B is the Boltzmann's constant, and T is the temperature in kelvin. The bulk ionic concentration n 0 can be obtained as where c is the concentration of the electrolyte in moles and N A is Avogadro constant. Therefore, the net charge density can be defined as Using hyperbolic functions, the net charge density can be expressed as The nonlinear Poisson-Boltzmann equation, obtained by substituting Eqs. (7) into (3), allows us to determine the internal potential distribution of the ions inside the fluid as Fluid flow phenomena due to EO in microchannels can be modelled by using the modified incompressible Navier-Stokes equations that can be written as follows.
• Continuity equation where ρ the fluid density and u i the fluid velocity components. By using the ideal gas law, the speed of sound c can be expressed as where p is the pressure. The continuity equation can now be rewritten as follows For incompressible flows, the speed of sound may be replaced with an artificial compressibility parameter, β , (Nithiarasu 2003).

• Momentum equation
where τ ij are the deviatoric stress components, related to velocity as in which µ is the dynamic viscosity. The last term of Eq. (12) takes into account the driving force of EOF due to the interaction between the EDL potential and the external electric field.

Dimensionless form of the governing equations
The dimensionless form of the governing equations for forced flow can be obtained through the following nondimensional scales.
• Laplace equation where κ is known as Debye length and corresponds to the EDL characteristic thickness (Patankar and Hu 1998) and L ref is a reference length that can be derived as follows: where A fluid is the area of the channel occupied by the fluid, L chan is the length of the channel, A chan is the area of the channel, A obstr is the area of the channel occupied by the obstructions. For the sake of clarity these areas are shown in Fig. 2. Depending on A fluid an effective porosity can be defined, as follows: By using the effective porosity the reference length can be rewritten as where W chan is the width of the channel. In plain channels L ref corresponds to the channel width (i.e. the distance between the two walls), consistent with the quantity commonly used in previous works concerning EOF. However in channels with obstructions, these reduce the area of the channel occupied by the fluid taking into account the variation of internal potential distribution due to the static charge of solid particles in addition to that of channel walls. In unobstructed flow, flow rate has been found to increase as the ratio between channel width and the Debye length, κ · W, is in the range of 1-100 and to be constant beyond this value (Rice and Whitehead 1965;Yao and Santiago 2003). When the channel width is much larger than the Debye length, the drag effect due to electro-kinetic forces of EDL decreases, and as a consequence EOF is reduced. Thus the use of charged flow obstructions to enhance EDL potential distribution can increase the range of effectiveness of EOF-driven systems.
Based on the above-defined quantities and scales, the non-dimensional form of the Navier-Stokes equations for these problems can be written as: • Continuity equation where the artificial compressibility parameter, β, Nithiarasu (2003) is locally calculated at each node as (Massarotti et al. 2006): The constant η is assumed to be 0.5, u conv and u diff are the convective and diffusive velocities, respectively.
• Momentum equation where The above set of non-dimensional equations has been solved by using fully explicit artificial compressibility-based CBS scheme (Nithiarasu 2003;Nithiarasu et al. 2016), described in "Appendix" for the sake of completeness.

Results and discussion
A silicon microchannel of 30 µm in width, characterised by an aspect ratio of 10, with deionised water as working fluid is considered as a reference case. The electrical conductivity, σ, and dynamic viscosity, µ, are assumed to be constant. The parameters used in the present study are reported in Table 1.
The flow past obstructions are investigated at the pore level, and the solid particles are assumed to be circular in shape. By fixing the effective porosity equal to 0.8, three different particle sizes are considered, i.e. particles with Fig. 2 Schematic of microchannel with obstructions diameter equal to 29, 16 and 12% of the channel width, hereafter, referred to as P29%, P16% and P12%, respectively. The layouts of these configurations of microchannel are shown in Fig. 3. The boundary conditions applied in this work are shown in Fig. 4.
The channel walls and the boundaries of solid obstructions are assumed to be active with a prescribed nondimensional zeta potential and to obey no-slip velocity boundary conditions. An applied external potential difference between the inlet and outlet is considered, and the normal components of velocity gradients are assumed to be zero at both inlet and outlet. The computation is started with prescribed zero velocity components as  initial condition. The dimensional external potential difference between the inlet and outlet sections of the channel is fixed equal to 1 kV/m. The non-dimensional zeta potential value of −0.75, corresponding to −19 mV, is considered for both the channel walls and the solid particles. The value considered for the zeta potential imposed on the charged surfaces is derived from some experimental investigations carried out on electro-osmotic flow in silicon microchannels (Eng 2009;Eng and Nithiarasu 2009). A set of 2D unstructured meshes, refined near all solid boundaries to capture the rapid change in both internal potential and velocity, is used. The details of the meshes used for plain channels and channels with obstructions are shown in Fig. 5. A mesh sensitivity study has been carried out to finalise the meshes used in the calculations.

Comparison of results against analytical solutions
The numerical model described in "Appendix" is used to determine EOF in microchannels with and without solid obstructions. For a plain channel without obstructions, the flow results are compared against the available analytical solutions (Patankar and Hu 1998). For a two-dimensional rectangular channel, the analytical solution for the internal potential may be written as where y is the distance from the wall. The analytical solution for the horizontal velocity component (Arnold 2007) may be written as In Fig. 6 the internal potential distribution and velocity profiles over half of the channel width are plotted and compared against the analytical solution. Since the analytical solution is derived for the linearised Poisson-Boltzmann equation, there is a slight discrepancy with the numerical results. The data of the simulations for velocity have been normalised in order to have a clear comparison with the analytical solution.

Effect of obstructions on EOF
The internal potential distribution plays a fundamental role on EOF. The influence is highlighted in Fig. 7, where the profiles of internal potential and velocity are plotted, over The internal potential and horizontal velocity profiles are constant along the channel length for plain channel, while they change when obstructions are introduced. This is shown by plotting the quantities of interest at different sections of the channel: in the central section (i.e. x/L = 5), where the profile presents a discontinuity due to the obstruction of a particle, and at the outlet section (i.e. x/L = 10) with no obstructing particle. Although the internal potential is equal to zero in a larger section for the plain channel, a larger internal potential gradient near the wall produces a higher horizontal velocity than that of the channel with obstructions. The solid obstructions enhance the internal potential distribution thanks to their charged boundaries, but at the same time they increase the resistance to fluid flow.

Effect of microchannel width on EOF
In order to assess the conditions under which the use of flow obstructions positively affects EOF, the geometry has been investigated. A range of different channel widths, W, between 5 and 150 µm, has been investigated for microchannels with and without obstructions. As the channel width increases, the internal potential approaches zero in the central region of the channel, and the profile becomes steeper close to the walls. This behaviour is weaker in channels filled with obstructions in which the presence of charged particles balances the width effect, as shown in Fig. 8. In this figure, the internal potential profiles at the outlet of the channel for a channel without and with obstructions are shown at different channel widths. Both internal potential and width are non-dimensional, and only half width of the channel is presented. The effects of obstructions on EOF can be appreciated in Fig. 9, where the non-dimensional horizontal velocity profile at the outlet section of the channel is reported for plain microchannel and microchannel with obstructions.
Since internal potential effect on EOF is introduced through a source term in momentum equation, confinement of internal potential variation close to walls produces little flow. For a plain fluid channel with no obstructions, the average velocity increases as the channel width is increased from 5 to 15µm, while beyond this value it rapidly decreases, as shown in Fig. 9b. This increase in velocity, initially, is a result of decrease in ratio of surface to cross-sectional areas, while EDL thickness remains constant. When flow obstructions are introduced in the microchannels, the average velocity increases beyond a channel width of 15 µm, indicating that the introduction of obstructions is effective beyond this channel width. However, the average velocity starts to decrease beyond a width of 60 µm, indicating that the interaction between the EDL of channel walls and obstruction boundaries becomes weaker in producing higher EO velocities. In addition, by comparing EOF in microchannels with and without obstructions, it is worth noticing that in plain channels the average velocity is higher than that of channels with obstructions for smaller microchannels, but it becomes significantly lower when the channel width is larger than 100 µm. Moreover, the rate of decrease in velocity in channels with obstructions is lower than that observed in plain channels.
As expected, in narrow channels, introducing obstructions produces higher drag resistance that undermines pressure generated by EO force. In addition EDLs overlap due to the short distance between adjacent particles. However, as the channel width is increased this trend is reversed and introduction of obstructions becomes more and more effective. On the other hand, when the width is too large, and the distance between the charged surfaces increases, the internal potential approaches zero in between obstructions, causing the average velocity to decrease. This effect is clearly noticed in Fig. 10, where the distribution of internal potential and horizontal velocity is reported for two different channel widths. Only part of the whole domain is plotted here to clearly view these details.

Width beyond which the solid obstruction is useful
In order to assess the effectiveness of using obstructions to enhance EOF, the flow rate has been determined for different channel widths and it is plotted in Fig. 11. From this figure it can be seen that for smaller plain channels the flow rate significantly grows as the width is increased. In all cases a rapid increase is observed at lower channel widths and the rate of increase is reduced beyond a channel width of 50 µm . The plain microchannel without obstructions appears to be effective up to a channel width of 80 µm, while beyond a width of 100 µm, introducing flow obstructions enhances the flow rate further. As shown in Fig. 11 three geometric configurations have been studied. By keeping the effective porosity constant, the ratio between the channel width and the particle diameter is varied and assumed equal to 29, 16 and 12%. Thus three different particle distributions within the category of channel with obstructions are analysed. They are identified, respectively, as P29%, P16% and P12%. All the configurations enhance EOF as the channel width is increased, but P16% appears to be more effective. The flow rate increases up to a channel width of 120 µm when this distribution is used, while smaller and larger particles enhance flow only up to around 90 and 100µm, respectively. The plain microchannels with no obstructions may be effective in pumping fluids up to a channel width of 100 µm. Beyond this channel width however, introduction of obstructions can increase fluid flow.
In the previous works concerning EOF, the velocity is often related to the EDL characteristic thickness, κ −1 . Under the analysed conditions this thickness is around 1µm. When the quantity W · κ approaches 100, the velocity is expected to rapidly decrease, as observed by several authors (Rice and Whitehead 1965;Yao and Santiago 2003). This is however not the case for channels with obstructions. The value of W · κ is much larger than 100 before a decrease in velocity is observed. If however L ref from Eq. 16 is used to define the width, W, for channels with obstructions, the flow rate decreases as W · κ approaches 100.

Effect of particle sizes on EOF
The effect of particle size on EOF appears to be significant. Configurations P29% and P16% present a similar trend in terms of flow rate for channel widths of up to 80 µm: beyond this value it appears that for P29% EDL effect on flow is reduced. This is due to the increased distance between the particles that becomes much higher than EDL thickness. For smaller particles (P12%) the flow rate is lower due to the lower distance between obstructions and EDL overlapping. The average velocity achievable in such microchannels is very small. This effect is shown in Fig. 12, where the EDL potential and horizontal velocity distributions are shown for a microchannel 60 µm in width. It appears that the decrease in particle to channel size ratio has a profound effect on flow rate. In particular a reduction in this quantity appears to reduce flow rate.

Effect of porosity on EOF
Finally, the influence of effective porosity on EOF is examined. Figure 13 shows the effect of porosity on flow enhancement when channel width is varied. The analysis is carried out in a high porosity range, from 0.4 to 0.8, since the results are finalised to be used for energy applications (Wang 2012). As seen a decrease in the porosity Flow rate at different widths of microchannels with and without obstructions from 0.8 to 0.4 increases the channel width range in which EO flow is enhanced. While flow has come to a stand still at a channel width of 240 µm when a porosity of 0.8 was employed, a finite flow rate is maintained at a porosity value of 0.6 and 0.4. It is also easy to notice that the flow pattern is approaching that of a channel with no obstructions when porosity is increased from 0.4 to 0.8. The effect of porosity largely confirms previous findings. As expected, in narrower channels the effect of overlapping EDL is more evident than in those characterised by a higher porosity, as highlighted in Fig. 14.

Concluding remarks
Electro-osmotic flow through microchannels with and without obstructions has been investigated. For modelling the flow through microchannels with obstructions, an equivalent reference length is introduced to represent the effect of obstructions. The model developed has been used to investigate the effectiveness of introducing flow obstructions in microchannels on enhancing electro-osmotic flow. The results show that beyond a channel width of 100 µm, introducing obstructions increases the effectiveness of electro-osmotic flow-driven systems. The results indicate that obstruction size affects the electro-osmotic flow and that an appropriate compromise between obstruction diameter and their distribution can be found in order to maximise flow rate. Finally, for larger microchannels, increasing the porosity decreases the flow rate. Thus, a decrease in porosity is recommended as the channel width is increased. The results of this study clearly show that the use of obstructions leads to flow enhancement in electro-osmotic systems for most cases analysed.
2. calculation of the pressure 3. correction of the velocity The resulting set of PDEs is then discretised in time by using a characteristic approach as it is described in the following subsections.

Temporal discretisation: CBS algorithm
In step 1, an intermediate velocity, ũ, is calculated after neglecting the pressure term.
At step 2 the pressure is determined as follows and it is used to correct the momentum in step 3 as Substituting u n+1 from (27) into Eq. (26), the final temporally discrete form of the pressure equation is obtained as where 0.5 < θ < 1.0 and θ 2 = 0 for a fully explicit scheme.

Spatial discretisation
The above equations are spatially discretised using the Galerkin finite element procedure and the following spatial discretisation of the variables.
where the overline indicates the nodal value, m is the number of nodes in the element, and a represents a specific node. After spatial discretisation the equations are weighted through the shape functions, N, and then integrated over the whole domain.
(26) The final matrix forms of the spatially discretised equations may be written as: (32) �u i = �ũ i + M −1 �t G p n + θ 2 �p − �t 2 Lp n L ≡ � ∂ ∂x i N T u i ∂N ∂x j d�;