Numerical study of dc-biased ac-electrokinetic flow over symmetrical electrodes

Numerical study of dc-biased ac-electrokinetic flow over symmetrical electrodes Wee Yang Ng, 1,2,a) Antonio Ramos, 3,b) Yee Cheong Lam, 2 and Isabel Rodriguez 1 1 Institute of Materials Research and Engineering, Agency for Science, Technology and Research (A*STAR), 3 Research Link, Singapore 117602, Singapore 2 School of Mechanical and Aerospace Engineering, Nanyang Technological University, Nanyang Avenue 50, Singapore 639798, Singapore 3 Departamento de Electro´nica y Electromagnetismo, Facultad de Fı´sica, Universidad de Sevilla, Avenida Reina Mercedes s/n, 41012 Sevilla, Spain (Received 6 August 2011; accepted 20 November 2011; published online 15 March 2012) This paper presents a numerical study of DC-biased AC-electrokinetic (DC-biased ACEK) flow over a pair of symmetrical electrodes. The flow mechanism is based on a transverse conductivity gradient created through incipient Faradaic reactions occurring at the electrodes when a DC-bias is applied. The DC biased AC electric field acting on this gradient generates a fluid flow in the form of vortexes. To understand more in depth the DC-biased ACEK flow mechanism, a phenomenological model is developed to study the effects of voltage, conductivity ratio, channel width, depth, and aspect ratio on the induced flow characteristics. It was found that flow velocity on the order of mm/s can be produced at higher voltage and conductivity ratio. Such rapid flow velocity is one of the highest reported in microsystems technology using electrokinetics. V C2012 American Institute of Physics. [doi:10.1063/1.3668262] I. INTRODUCTION Dynamic manipulation of fluids in microchannel is crucial for the development of lab-on-achip (LOC) devices. 1,2 One key requirement for the actuation system in LOC devices is to avoid mechanical moving parts which can potentially complicate the fabrication and increase the size of the system. 3,4 Electrokinetics has been demonstrated as a viable approach to effectively manipulate fluids in microchannels through the use of an electrical driving force. 5,6 One key advantage lies in its ease of implementation which makes electroosmosis (i.e., the use of DC electric field to drive fluid flow, due to the Coulomb force on the diffuse double layer) one of the most widely used pumping approaches in microfluidics in addition to pressure-driven methods. 7,8 AC-electrokinetic (ACEK) makes use of AC electric fields acting on AC charge density to induce rectified flow; two important examples of ACEK are AC-electroosmotic (ACEO) and AC-electrothermal (ACET) flows. 9–16 ACEO flow is based on the Coulomb force acting on the induced charge in the double layer in the presence of a tangential AC electric field. 10–13 ACET flow is based on the interaction of AC electric fields with conductivity gradients in the fluid induced by a thermal gradient. 14–16 More recently, a DC-bias AC voltage (V applied ¼V DC þV AC cos xt) has been employed effectively to concentrate particles/cells 17–19 and for pumping/mixing 20–23 applications in microfluidics. There has been several hypothesis on the driving mechanism of the flow induced by DC-bias AC voltage. 6,18,19 Recently, we conducted experimental investigation to elucidate the driving mechanism of the DC-biased ACEK flow phenomenon. 21,24 In the publications, we a) Electronic mail: [email protected]. b) Electronic mail: [email protected]. 1932-1058/2012/6(1)/012817/10/$30.00 V C2012 American Institute of Physics6, 012817-1 BIOMICROFLUIDICS 6, 012817 (2012) demonstrated that the applied DC-bias voltage gives rise to incipient Faradaic reactions at the electrodes’ surface. 21,24 Indeed, it was tested that these reactions cause changes on the local pH due to the generation of co-ions above the electrodes. As a result, the region above the anode became acidic due to the increase of H þ ions, and the region above the cathode was basic due to the increase of OH – ions. In addition, it was measured that the basic region is 11%–12% more conductive than the acidic region. 24 Hitherto, there is no existing theoretical model to describe the DC-biased ACEK flow. In this paper, we present a phenomenological model based on the actual conductivity difference between the anode and the cathode to give a description of the fluid flow. It should be noted that the effect of DC-bias could not be included by simply adding a DC-bias to the AC voltage in established ACEO model. 12,13 This is because we do not know with certainty all the electrochemical reactions taking place at the electrodes’ surface. Therefore, we take the value of the increment of conductivity from experiments. 24 Moreover, since the DC voltage is just enough to generate the species, the applied DC voltage is dropped across the double layer and not in the bulk. A convection-diffusion equation is used here to describe the distribution of conductivity. The interaction of the electric field with the conductivity gradient leads to fluid flow. This paper will demonstrate that the developed model could predict well the experimental observations 23,24 namely the single and dominant flow vortex and the direction of vortex rotation. It is experimentally difficult to characterize the cross-sectional fluid velocity profiles; therefore a numerical simulation serves as a tool to better understand the various factors affecting the induced flow velocity. Parametric factors (voltage, conductivity ratio) and geometrical factors (width, depth, and aspect ratio) were varied in this study to understand their effects on the induced DC-biased ACEK flow velocity. II. NUMERICAL ANALYSIS DC-biased ACEK flow phenomenon was modeled using COMSOL MULTIPHYSICS V R . A phenomenological model was presented here to give a description of the DC-biased AC-electrokinetic fluid flow. The numerical model is based on the measured conductivity difference between the anode and the cathode and was used as the basis for the boundary conditions for the conductivity equation. It was determined that the liquid above the cathode is more conductive than the one on the anode side. AC impedance analysis was adopted to measure the conductivity level above the cathode and the anode. It was found that with a DC bias of 2 V DC , the cathode is more conductive than the anode by 11%–12%. 24 We assume that the applied DC voltage is just above the threshold to generate Faradaic reactions at the electrodes. In this way, V DC is dropped at the electrode/electrolyte interface to drive the reactions and negligible DC voltage is dropped across the bulk electrolyte. Notice that the minimum voltage to produce electrolysis of water is about 1.23 V. Therefore, we can neglect V DC in the medium and solve only for the AC potential (V applied ¼V AC cos xt). Moreover, since the frequency of the AC signal is much greater than the reciprocal RC time for charging the double layer, 10 the AC voltage is dropped mainly across the bulk electrolyte and negligible AC voltage is dropped at the electrode/electrolyte double layers. In addition, the effect of joule heating is ignored since the conductivity used is low. In this study, a 2D analysis was performed for a pair of coplanar and symmetrical electrodes enclosed by a microchannel as shown in Fig. 1. To limit the scope in this study, the electrode width and gap were fixed at 40 lm and 20 lm, respectively, but the microchannel width (W) and depth (D) were varied from 100 to 400 lm and 20 to 200 lm, respectively. A. Electrical equations For an electrolyte, the electric potential (/) distribution at any point in the solution can be related to the electric charge density (qe) which is described by the Poisson equation rðe~ EÞ¼qe;(1) 012817-2 Ng et al. Biomicrofluidics 6, 012817 (2012) where eis the permittivity and ~ E¼r/is the electric field. The movement of the ions in the solution will constitute to an electrical current ~ J¼Xqinili~ EqiDirniþqini~ V  ;(2) where qi(¼zie, with zithe ionic valence and ethe elementary charge) is the ion charge, niis the ion number density, li(¼zieD i=kT) is the mobility, Diis the diffusion coefficient, and ~ Vis the velocity. The first, second, and third terms represent, respectively, electromigration (i.e., conduction), diffusion, and convection. Typically for electrolytes, the convective current qini~ V is usually smaller than the conduction current qiniliEwhich is governed by the electrical Reynolds number. 25 In microsystems, the electrical Reynolds number is very small for electrolyte. For semi-insulating liquids, a conductivity smaller than 10 9 S/m would be required in order to have a convective current of the same order than the conduction current. In addition, the ratio of diffusive current to electromigration Dirni jj=jnili~ Ejis smaller than unity in the liquid bulk. 25 Typically, for voltage drop in the bulk much greater than 0.025 V, diffusion is negligible as compared to electromigration. The ratio between these two currents is on the order of 0.025/V. In our experiments, we apply around 10 V 0.025 V. In this situation, we can write ~ J¼r~ E, where r¼Pqiniliis the electrolyte conductivity. With these approximations, the charge conservation equation for the AC potential becomes 25–27 rððrþiexÞrUÞ¼0;(3) where Uis the phasor of the electrical potential. Consider a binary electrolyte (ni¼1), assuming electroneutrality, 16 the solution conductivity can be represented in the model with a convection–diffusion equation 27–29 @r @tDr2rþð~ VrÞr¼0;(4) where the second and third terms represent, respectively, diffusion and convection. B. Mechanical equations In a microscale system, Reynolds number is usually small (<1) and flow is laminar in nature. The fluid velocity is governed by the Navier-Stokes equations, for steady-state solutions, we have qfð~ VrÞ~ V¼rPþgr2~ Vþ~ Fe;(5) FIG. 1. Geometrical configuration and boundary conditions of the simulated DC-biased AC-electrokinetics problem. Boundary conditions for right electrode (conductivity: r Right , electrical: V GND , flow: no slip), left electrode (conductivity: r Left , electrical: V AC , flow: no slip), and channel walls (conductivity: insulation=symmetry, electrical: electrical insulation flow: no slip). 012817-3 dc-biased ac-electrokinetic flow Biomicrofluidics 6, 012817 (2012) r~ V¼0;(6) where qfis the fluid density, Pis the pressure, gis the dynamic viscosity, and ~ Feis the electrical body force. The electrical body force on the fluid is given by the time-average Coulomb force 25–27 ~ Fe¼qe~ E  ¼1 2Re½qe~ E;(7) where qeis the charge density phasor, Re[...] represents the real part, ~ Eis the complex conjugate of the electric field phasor. In this case, the time-average Coulomb force is used in the calculation. This is because the ac signals have periods much shorter than the typical mechanical time of variation of the system. Despite the electrolyte being quasi-electroneutral, the residual charge can lead to a significant force. The charge density (qe) in the bulk can be obtained from Eqs. (1) and (3) leading to a charge density phasor 26,27 qe¼e rþiex rrrU:(8) C. Numerical method and settings To solve for the DC-biased ACEK flow, “convection and diffusion,” “conductive media DC,” and “incompressible Navier-Stokes” modules in the COMSOL software were employed. COMSOL software employs the finite element method (FEM) in solving the partial differential equations (PDEs). The numerical scheme discretizes the defined spatial domain into triangular finite elements to approximate the solutions of the PDEs. Typical computation consists of approximately 1300 triangular elements and the minimum element quality is 0.86. Further investigation on the numerical convergence was conducted. With mesh refinement up to 5000 or even 20 000 triangular elements, the computed velocity does not change appreciably (less than 1%). In the module of “convection and diffusion,” constant conductivities (r1and r2) were set on the two electrodes and insulation/symmetry on the channel walls. In the module of “conductive media DC,” various potentials were set on the two electrodes and electrical insulation was applied on the channel walls. In the module of “incompressible Navier-Stokes,” no slip boundary conditions were set for the electrodes and the channel walls. In the solving process, an iterative approach was adopted where the problem was solved in three steps: First, the convection-diffusion equation was solved for the conductivity field, Eq. (4). Next, the computed conductivity solution was used to solve for the electrical potential field, Eq. (3). A parametric solver was adopted to set the voltage where the voltage was increased from zero to the set voltage level at a voltage step of 0.01 V. Lastly, the computed potential solution was used to solve for the electrokinetic flow velocity field, Eqs. (5) and (6). The solving process was repeated until we obtained a converged solution where the relative tolerance is set at 10 3 . The parameters used in the simulation are listed in Table I. III. RESULTS AND DISCUSSIONS A. DC-biased AC-electrokinetic flow Figure 2shows the simulated results for the velocity distribution in the 2D cross-section model (W ¼100 lm, D ¼50 lm, electrode width ¼40 lm, and electrode gap ¼20 lm) for both negative and positive DC-biased AC-electrokinetic flow. The voltage amplitude used was 10 V and the conductivity level at the two electrodes was set with a 10% difference. In both cases, a single and dominant flow vortex was obtained in the direction from higher conductivity (basic, cathode) to lower conductivity (acidic, anode). As explained, due to Faradaic reactions, a conductivity gradient is generated from the basic (higher conductivity, net negative electrode) to 012817-4 Ng et al. Biomicrofluidics 6, 012817 (2012) the acidic (lower conductivity, net positive electrode) domains. When an AC signal is applied, a horizontal force component is generated close to the electrode surface where the field is strongest. As such, the fluid is actuated from the higher conductivity to the lower conductivity domains close to the surface of the electrodes; this determines the direction of rotation and eventually results in a single and dominant vortex over the electrodes pair. 24 For a negative DC-bias AC voltage case (see Fig. 2(a)), the fluid flows from right to left, i.e., from higher conductivity OH – ions rich cathode to lower conductivity H þ ions rich anode, at region near the electrodes’ surface. This flow produces a resultant clockwise (CW) vortex rotation. For a positive DC-bias AC voltage case (see Fig. 2(b)), the fluid flows from left to right, i.e., it is still flowing from higher conductivity OH – ions rich cathode to lower conductivity H þ ions rich anode, at region near the electrodes’ surfaces. This flow produces a resultant anti-clockwise (anti-CW) vortex rotation. As such, the predicted results on the vortex direction corroborate the reported experimental observations. 23,24 TABLE I. Simulation parameters. Property Value Fluid density (q)110 3 kg=m 3 Fluid dynamic viscosity (g)110 3 kg=ms Relative permittivity of KCl (er)80 Permittivity of free space (e0) 8.85 10 12 F=m Coefficient of diffusion (D)110 9 m 2 =s Voltage amplitude (V) 1–20 V AC frequency (f) 100 kHz Conductivity at electrode 1 (r1) 1–10 mS=m Conductivity at electrode 2 (r2) 1–10 mS=m FIG. 2. Simulated velocity distributions. (a) Negative DC-bias, r1<r2and (b) positive DC-bias, r1>r2. Conditions: voltage amplitude ¼10 V, conductivity difference ¼10%, W ¼100 lm, D ¼50 lm, electrode width ¼40 lm, and electrode gap ¼20 lm. For negative DC-bias of (a), plots of X-velocity (horizontal) across (c) channel depth at various locations along the left electrode, and (d) channel width at different channel heights. 012817-5 dc-biased ac-electrokinetic flow Biomicrofluidics 6, 012817 (2012) In general, the negative and the positive DC-bias simulated results are similar and symmetrical. The x-velocities (Vx) of the simulated results for the negative DC-biased case (Fig. 2(a)) were extracted and plotted. The horizontal velocities were plotted across the width and across the depth at different locations on the channel as shown in Figs. 2(c) and 2(d), respectively. From Figs. 2(c) and 2(d), it can be observed that Vxare faster on the region nearer to the electrode surface. This is expected as the velocity is driven by the gradients of conductivity and the electrical potential, where the conductivity and the electric field are the strongest near the surface at the electrode gap. Taking the point where Vx¼0 along the center of the channel (y ¼0) as the center of the vortex; the center of vortex is at a height of 20 lm. This corresponds to 40% (i.e., 20 lm/ 50 lm) of the total channel height. Defining H v as the percentage of the vortex height to the channel depth (D), the H v value computed is thus rather high at 40%. Similarly, we have previously reported an observed high H v value (>30%) by experimental imaging the induced DC-biased ACEK flow vortex. 24 This high center of vortex implies that the flow mechanism is bulk driven. In the literature, ACET flow which is bulk induced generally has a vortex center between 30% and 40% of the channel height. 14–16 B. Effect of voltage and conductivity ratio The applied voltage and the conductivity ratio affect the velocities of the flow induced in the system. Here, the conductivity ratio (c) is defined as the ratio of the higher conductivity (r2) to the lower conductivity (r1). To quantify the velocities of the flow at various conditions (voltage, conductivity ratio), the average velocity amplitude in the cross-section domain was computed and defined as 27 Vave ¼1 AððVdxdy;(9) where V¼jVðx;yÞj is the velocity amplitude at ðx;yÞand A(m 2 ) is the area of the crosssection domain. Similarly, the adopted geometry to be study were W ¼100 lm, D ¼50 lm, electrode width ¼40 lm, and electrode gap ¼20 lm. The average velocity amplitudes for negative DC-bias at different parameters (voltage, conductivity ratio) were calculated and presented in Fig. 3. Figure 3(a) shows that the average velocity amplitudes increase with the applied voltage at various conductivity ratios. The highest velocities obtained for c¼1.1, 1.3, and 1.5 at 20 V were 2.78 mm/s, 7.64 mm/s, and 11.78 mm/s, respectively. Next, by fixing the voltage amplitude at 10 V, the effect of conductivity ratio on velocity was investigated. Figure 3depicts results for various conductivity ratios. It indicates a clear trend that a higher applied voltage and/or higher conductivity ratio will result in a faster flow velocity. From Fig. 3(a), the simulated velocity at 10 V and c¼1.1 is 695 lm/s, which compares favorably with the experimental measured velocity of 700 lm/s. 24 It should be noted that the conductivity gradient/ratio generated by the DC-bias Faradaic charging has a limitation. This is because electrolysis and electrodes degradation will occur at too high voltage. Therefore, the conductivity ratio via this method is limited by the level of DC-bias (i.e., <3V DC ). Nevertheless, to facilitate an understanding of mechanism, simulation is conducted for conductivity ratio up to 10 which is the typical ratio adopted in electrokinetic instability (EKI) flow phenomenon. 30–32 From Fig. 3(d), the achieved velocity at 20 V and c¼10 was extremely fast, 14.44 mm/s. As a comparison with a more complex numerical model, 33 it was reported that the achieved velocity was up to tens of mm/s (i.e., at 0.5 10 5 V/m and c¼10). Their predicted velocity is of the same order of magnitude as obtained by the current model. In addition, this result implies that rapid induced flow can be achieved by incorporating a microelectrode (with applied AC voltage) with manually induced conductivity streams. Typically, the EKI configuration introduces the high AC voltage through 012817-6 Ng et al. Biomicrofluidics 6, 012817 (2012) bulk electrodes from the reservoirs. Integrated microelectrodes can indeed supply electric field of the same order of magnitude (10 4 –10 5 V/m) due to the reduced electrodes gap. C. Effect of channel width, depth, and aspect ratio Channel dimensions can affect the effectiveness of fluid manipulation in the microchannel. As such, the effects of the channel width (W), depth (D), and aspect ratio (i.e., AR ¼D/W) on the average vortex flow velocity (Vave) were investigated. Similarly, the electrode width and gap were fixed at 40 lm and 20 lm, respectively. The voltage amplitude used was 10 V and conductivity ratio was at c¼1.1. Figure 4(a) shows the average velocity amplitudes (Vave)as the width of the channel (W) varies from 100 lmto400lm, at fixed depth of either 50 lmor 100 lm. It can be observed that the Vave generally decreases as the W increases. This is because there is a limit where the induced DC-biased ACEK flow can drag the surrounding fluids into motion. 34 Therefore, there is an optimum range where the fluids can be effectively actuated, determined by the driving velocity (i.e., velocity near the surface of the electrodes). As W increases (see Fig. 4(a)), the velocities are only confined within the central region of the microchannel, and thus, the actuation efficient generally decreases. These findings have implication on the design of pump and mixer devices, as it would be preferred to locate the electrodes in close proximity of the pumping/mixing region. The inset in Fig. 4(a) shows the corresponding H v percentage of the simulated data points. The trend indicates that the vortex center height stays relatively constant (28%–40%) for different channel width (100–400 lm) investigated. This is because of the fixed channel depths (50 lm, 100 lm) adopted. Figure 4(b) shows the average velocity amplitudes (Vave) as the depth of the channel (D) varies from 20 lmto200lm, at a fixed width of either 100 lm or 200 lm. The average velocity increases to a maximum and then decreases with increasing channel depth. It can be observed that there is an optimal depth where a maximum average velocity is obtained. The maximum velocities were 748 lm/s (at W ¼100 lm, D ¼100 lm, AR ¼1) and 698 lm/s (at FIG. 3. Average velocity amplitude at different (a) applied voltage, (b) conductivity ratio, up to 1.1, (c) conductivity ratio, up to 2, and (d) conductivity ratio, up to 10. Configurations: W ¼100 lm, D ¼50 lm, electrode width ¼40 lm, and electrode gap ¼20 lm. 012817-7 dc-biased ac-electrokinetic flow Biomicrofluidics 6, 012817 (2012) W¼200 lm, D ¼150 lm, AR ¼3/4) for the fixed widths investigated. Note that when the cross-section is large, the fluid velocity does not extend very far from the gap between electrodes. This means that the average velocity should decrease, since regions far from the gap has very small velocity. If the cross-section is small (short height), the velocity is reduced because of viscous friction with the walls. Therefore, there is a maximum in between the two limits. In addition, the maximum average velocity generally occurs at channel AR close to 1, this gives us clue on designing the geometry of the pump with better efficiency. The result indicates that for large channel cross-section (AR 1) the induced flow would not extend far from the electrode surface, hence will not be efficient. For small channel cross-section (AR 1), the induced flow velocity is decreased from expected due to viscous friction. Next, the inset in Fig. 4(a) shows the corresponding H v percentage of the simulated data points. The H v percentage generally decreases with increasing depth. Similarly, this is because the induced DC-biased ACEK flow also has a limit in setting the fluids in motion in the ydirection (depth-wise). In Fig. 4(b), it should be noted that for the fixed width of 200 lm, the channel AR <1, and the obtained H v percentage is still considerably high (>23%). For the fixed width of 100 lm, as the geometrical configuration changes from a flat microchannel (AR <1) to a tall microchannel (AR >1), a dramatic drop in the H v percentage is observed. This also suggests that the design of microchannel geometry should adopt a square (AR ¼1) or flat (AR <1) type of channel cross-section for electrokinetic flow actuation to be more effective. FIG. 4. Average velocity amplitude with varying channel (a) width, for a fixed depth of 50 lm or 100 lm, and (b) depth, for a fixed width of 100 lm or 200 lm. Inset shows H v value of each data point. Conditions: voltage amplitude ¼10 V, conductivity difference ¼10% (c¼1.1), electrode width ¼40 lm, and electrode gap ¼20 lm. 012817-8 Ng et al. Biomicrofluidics 6, 012817 (2012) As observed in Fig. 4(b), the channel AR can affect the manipulation of fluid inside a microchannel. Here, we investigate the velocity of the induced flow in a channel of various ARs (1/2, 1/3, 1/4, and 1/5), but we restricted our study to AR <1 (flat microchannel), the preferred channel configuration. Figure 5shows the average velocity amplitudes in different channel widths for various constant channel ARs. Figs. 5(a)–5(d) show the plots for AR of 1/2, 1/3, 1/4, and 1/5, respectively. With a decrease in AR from 1/2 to 1/5, the Vave generally decreases. The AR measures the flatness of the microchannel, thus a lower AR indicates that W is much larger than the D. Note in Fig. 4(a) with fixed D, the computed average velocity is also found to decrease with increasing W. Therefore, the trend in Figs. 4(a) and 5are the same depicting a decreasing AR with fixed D and increasing W. Typically in such configuration, the induced DC-biased ACEK flow has a limited actuation range as discussed previously for Fig. 4(a). It should be noted that the H v percentage remains relatively high at 27%–39%, 30%–42%, 33%–44%, and 35%–45% for AR of 1/2, 1/3, 1/4, and 1/5. IV. CONCLUSIONS In this paper, a numerical model describing DC-biased AC-electrokinetic flow has been presented. “Convection and diffusion,” “conductive media DC,” and “incompressible NavierStokes” modules were adopted in the COMSOL software. The phenomenological model is based on experimentally measured conductivity gradient generated upon the application of a DC-bias to one of the electrode pairs. The simulated results indicate good correlation with reported experiments in terms of flow induced vortex characteristics, magnitude, and rotating direction. It was confirmed that the flow direction was from the higher conductivity to the lower conductivity regions near the surfaces of the electrodes and thus forming a single and unified vortex across the width of the microchannel. Parametric studies were conducted by varying the voltage, conductivity ratio, channel width, depth and aspect ratio. It was found that the flow velocity induced in DC-ACEK can reach values in the range mm/s by selecting an appropriated higher voltage and conductivity ratio. In addition, it was found that micro channels with square dimensions (AR ¼1) are more efficient for optimum flow velocity. FIG. 5. Average velocity amplitude at different channel aspect ratios (a) AR ¼1=2, (b) AR ¼1=3, (c) AR ¼1=4, and (d) AR ¼1=5. Inset shows H v value of each data point. Conditions: voltage amplitude ¼10 V, conductivity difference ¼10% (c¼1.1), electrode width ¼40 lm, and electrode gap ¼20 lm. 012817-9 dc-biased ac-electrokinetic flow Biomicrofluidics 6, 012817 (2012)