Electroosmotic flow: From microfluidics to nanofluidics

Abstract Electroosmotic flow (EOF), a consequence of an imposed electric field onto an electrolyte solution in the tangential direction of a charged surface, has emerged as an important phenomenon in electrokinetic transport at the micro/nanoscale. Because of their ability to efficiently pump liquids in miniaturized systems without incorporating any mechanical parts, electroosmotic methods for fluid pumping have been adopted in versatile applications—from biotechnology to environmental science. To understand the electrokinetic pumping mechanism, it is crucial to identify the role of an ionically polarized layer, the so‐called electrical double layer (EDL), which forms in the vicinity of a charged solid–liquid interface, as well as the characteristic length scale of the conducting media. Therefore, in this tutorial review, we summarize the development of electrical double layer models from a historical point of view to elucidate the interplay and configuration of water molecules and ions in the vicinity of a solid–liquid interface. Moreover, we discuss the physicochemical phenomena owing to the interaction of electrical double layer when the characteristic length of the conducting media is decreased from the microscale to the nanoscale. Finally, we highlight the pioneering studies and the most recent works on electro osmotic flow devoted to both theoretical and experimental aspects.


Introduction
Ever since Ferdinand Friedrich Reuss reported for the first time, more than two centuries ago, his interesting observation of moving water through a plug of clay upon application of an external electric field, a huge body of theoretical and experimental studies have been carried out based upon the discovery. The most important outcome of this discovery was the ability to make water flow without any mechanical parts and solely through application of an external electric field. The working principle of this concept, which was later called electrokinetic transport, involves the interplay of the external electric field and a charged interface that is neutralized by a counter-charge layer of the liquid. One of the well-known categories of electrokinetic transport is EOF. EOF, which can be employed as a pump to drive water from one side to the other side of an electrically charged medium, has been extensively used and studied [1][2][3][4][5][6][7][8][9][10][11]. The key role in this pumping phenomenon is played by the electrically charged solid-liquid interface, the so-called electrical double layer (EDL). The EDL, which is a composition of charged solid surface and a very tiny layer of counter-charges in an aqueous solution (around a few nanometers), was first investigated in the 19th century [12], and has been further studied in the 20th century [13][14][15] and in contemporary times [16][17][18][19][20]. Consequently, it is of utmost importance to enhance our understanding of the configuration of water and ionic species in the vicinity of an electrically charged solid surface.
One important characteristic feature of EOF is that its volume flow rate can be comparable to that of pressure-driven flow in micro-and nanoscale media [21]. Therefore, it is important to study EOF in micro-and nanochannels/porous media and understand the impact on the flow of water or charged species. It is worth pointing out that the characteristic size of the medium has a considerable impact on EOF due to the overlapping of EDLs. Let us imagine that the height of a channel is comparable with the EDL thickness. Thus, it is expected that the EDLs interact. This interaction of the EDLs Figure 1. Reproduced schematic diagram [28] of the experiments done by Reuss from his article [28]. The experiments conducted in a U-type glass tube filled with water and the lower part was filled with insoluble particles, such as sandstone, which created a porous barrier. Reuss observed that when an external voltage was applied to the water, it began to pass through the porous barrier from the anode (+) to the cathode (−) side. makes a major part of the channel electrically non-neutral [22,23]. As a result, due to the electrostatic forces, the channel will become selective to the counter-ions, repelling most of the coions. This interesting phenomenon has been the basis of many practical applications-from biophysics [24] to electrokinetic remediation (EKR) of contaminants from underground water resources [25,26].
This tutorial aims at reviewing EOF at the microfluidic and nanofluidic scale. We start from a short historical review of EOF and the story of its discovery (Section 2). After discussing the background, we present and discuss the evolution of the various EDL models (Section 2.1). After understanding how a solution in the vicinity of a charged solid surface will be configured and a chemically active solid surface will be electrically charged, we will then discuss EOF through microchannels (Section 3.1) and microporous media (Section 3.2). Next, we will turn our interest down to nanoscale channels (Section 4.1) and nanoporous media (Section 4.2), which we will address both theoretically and experimentally. We present out conclusions in .

History and background
When German scientist Ferdinand Friedrich Reuss reported his observation of water flowing through a plug of clay under application of an external electrical voltage through the two ends of a U-type tube, he presented his discovery as a lecture entitled "Notice on a new, hitherto unknown effect of galvanic electricity" to the Physico-Medical Society of Moscow in 1807 [27,28]. Subsequently, he published two papers in Latin and French [29,30] to describe this unknown phenomenon in detail. In his first experiment (Fig. 1), he observed the flow of water (electroosmosis) through a wetted porous barrier made of powdered quartz. Upon application of an external electric field, the water on the cathode (−) side was raised while that on the anode (+) side descended. After 2 h, he observed that  [28,30]. Reuss inserted two glass tubes into a block of moist clay (part A) and applied an external electric field to the water inside the glass tubes.
the entire cathode side of the U-tube, together with the S-type glass tube, was filled with water while the anode side of the U-tube was discharged.
Reuss conducted his second experiment [30] in which two glass tubes were filled with water and a thin layer of wellwashed sand was kept at the bottom (Fig. 2). In this setup, Reuss observed that by applying an external electric field, the sand at the bottom of the positive tube swelled upward and tended to penetrate the sand layer (demonstrated as label A in Fig. 2). With the benefit of hindsight, we today will interpret this phenomenon as electrophoretic movement of the charged sand particles. However, Reuss interpreted this observation as EOF into the sand.
Reuss made several false deductions while trying to interpret this phenomenon. He stated that the liquid between the poles of a battery is continuously driven from the positive pole toward the negative pole. He also believed that the presence of porous barriers makes this movement visible by counteracting the impact of gravity [28]. While some previous works supported his conclusions [31,32], later in this article we will show that according to our modern understanding, these conclusions are false. Biscombe [28] has reviewed the history of EOF discovery, which can be referred for further details.

Models to explain the charged solid-liquid interface
As we mentioned in the last section, Reuss concluded that the fluid flow must be due to the application of an external electric field in an aqueous solution. However, the main physics underlying this phenomenon was a mystery to him. He was not aware of the impact of the solid-liquid interface on it. In 1859, Georg Quincke [33] conducted several important experiments that proved to be a great step forward in shedding light on Reuss' discovery [32]. He conducted the reverse electroosmosis experiment by applying a pressure gradient. In his setup, the pumped water went through a tube where he measured a potential difference. He found that adding sodium chloride lowered the measured electric potential. The main conclusion from his experiments was that the sign of the electric potential is independent of the water flow rate, tube diameter, the concentration of the dissolved ions, and even the porous barrier material in which he employed the glass, sand, graphite, silk, and ivory. Here, we should note that, in his experiments, the measured electric potential was essentially changed by utilizing different materials. Quincke's experiments led him to the idea that there must be an excess space charge rather than the surface charge. This hypothesis was a remarkable milestone for the following works that attempted to explain the underlying physics of the earlier experiments. The pioneering theoretical attempt was done by the German physicist and physician Herman von Helmholtz in 1850 [12]. Helmholtz was the first who paid attention to the nature of the electrode-electrolyte or charged solid-electrolyte interface. Essentially, a charged interface will attract free counterions in the solution and repel coions because of the Columbic forces. The charged interface, along the layer of counterions that are attracted to the interface, was called the EDL. Helmholtz assumed that there would be no electron transfer at the charged interface and that the ions were solid spherical particles with a specific diameter. According to the Helmholtz model, as illustrated in Fig. 3, the charge that holds (electrode) or is acquired on the solid-liquid interface must be balanced by the redistribution of ions arranged in a parallel hypothetical plane with the solid-liquid interface. The development of the kinetic theory of molecular behavior proved that the Helmholtz EDL model was not realistic [34]. The main shortcoming of the Helmholtz model was in ignoring the thermal motion of the counter-ions and the possibility of adsorption onto the solid surface. In the Helmholtz model, the thickness of the Helmholtz double layer is independent of the thermochemical properties of the solution.
Taking into account the thermal motion of the ionic species, it is reasonable to say that the attracted counter-ions in the vicinity of the solid-liquid interface will spread out in space [21]. The spreading of the counter-ions into the bulk solution forms a diffuse layer that, contrary to the Helmholtz double layer, varies with the thermochemical properties of the solution. The idea of such a diffuse layer was independently proposed by French physicist Louis Georges Gouy (1910) and English physical chemist David Leonard Chapman (1913). The resulting model is now named after them as the Gouy-Chapman electrical double layer (GC-EDL) model. According to GC-EDL (Fig. 4) model, the electric potential decreases exponentially with the distance from the charged solid-liquid interface. The thickness of the region with nonzero net electric charge density ( ρ e = e(c + − c − ) where e denotes the electron charge and c + and c − are the counter-ion and co-ion concentration, respectively), depends on the bulk ionic concentration that, for instance, in the solution would be stretched over 100 nm in a very dilute electrolyte solution.
The standard solution of the GC-EDL model, which is based on statistical mechanics, can be found in textbooks as early as Adamson [35] and Overbeek [36]. The GC-EDL model is based on the Poisson equation, which is solved from the solid surface (at x = 0, ψ = −ζ ) to the bulk solution (where x → ∞ and ψ = 0) and is written as where ψ, ε 0 , ε r denote the electric potential, vacuum electrical permittivity, and the solution's relative permittivity with respect to vacuum, respectively. For the ionic species in the solution, it is reasonable to say that they feel the local electrostatic potential from the charged solid-liquid interface. Consequently, we can introduce the Boltzmann equation as where c i , c b i , z i , k B , and T represent the ith ionic concentration, the bulk ion concentration, the ionic valence (e.g., for a monovalent ionic solution z 1 = +1 and z 2 = −1), Boltzmann constant ( k B = 1.380649 × 10 −23 J/K), and the absolute temperature of the solution, respectively. Introducing Eq. (2) into Eq. (1) gives rise to the Poisson-Boltzmann equation: The right-hand side of Eq. (3) can be simplified if we assume that ψ k B T e . The term V t = k B T/e is well-known as the thermal potential and provides a measure of the induced potential energy on an elementary charge (i.e., electron charge) [37]. At room temperature, V t is about 25 mV. Recalling the assumption for the local electric potential in solution, if ψ 25 mV, then Eq. (3) can be simplified as which linearizes the Poisson-Boltzmann equation. Historically, Debye and Hückel simply extended the exponential term related to the Boltzmann equation as truncated Taylor series to the first order [38]. In the Taylor series, the zerothorder vanishes because the whole system is electroneutral ( n i z i ec i = 0). However, the first order leaves the Helmholtztype equation. Rewriting Eq. (4) by only keeping the Laplace term at the left-hand side, we finally obtain Considering the right-hand side of Eq. (5), the term κ = ( 2e 2 c b i ε 0 εr k B T ) has got the inverse of a length, which is attributed to the characteristic length of the EDL and called the Debye length. Consequently, Eq. (5) can be rewritten as It is worth pointing out that κ is also referred as the Debye-Hückel parameter. Solving Eq. (6) together with the Boltzmann equation (Eq. (2)) determines the distribution of the ionic species and the electric potential in the vicinity of a charged solid-liquid interface. Regarding the solution for Eq. (6), we can subject the equation to the following boundary conditions: a) At y = H, dψ dy = 0 due to symmetry at center of channel , b) At y = 0, ψ = ζ (potential at shear plane) , Figure 5. Reproduced sketch of the Stern EDL model [41] where the solid surface is assumed to be positively charged. In his model, the first layer, called Stern layer, has the thickness represented by δ, and the diffuse layer assumed to be beyond the Stern layer.
In Eq. (8), H represents the distance far from a single solid-liquid interface or in the middle of a 2D micro-or nanochannel. Here, it is worth pointing that we can also propose an analytical solution for the Poisson-Boltzmann equation (Gouy-Chapman theory, Eq. (3)) without any approximation (Debye-Hückel). It should be noted that the GC-EDL model is a nonlinear theory for ionic distribution in the vicinity of the solid-liquid interface. In this regard, the analytical solution for Eq.
where denotes the nondimensional form of ψ defined as = ψ Vt ; at x = 0 and x → ∞ we have = s and = 0, respectively. There are several excellent textbooks that can be referred for further study about the electrostatic interaction between planar and curved surfaces [21,39].
Following the Helmholtz and GC-EDL models, researchers attempted to propose more complicated EDL models that could tackle the shortcomings of these previous models. For instance, the Gouy-Chapman model suffers from a lack of generality for highly charged solid-liquid interfaces. According to the GC-EDL model, if we increase zeta potential to a very high value, the distributed ionic concentration is predicted to be infinitely large. This prediction originates from assuming the ionic species as point charges (no ionic volume effect). However, in reality, the ionic species have nonzero volume that approaches the solid surface to a distance not less than their Stokes radii or the effective hydrated radius in solution [40]. Thus, the German-American physicist Otto Stern (1924) proposed a combination of the Helmholtz and GC-EDL models to overcome the shortcomings of both models [41]. He suggested that the EDL in the solution consists of two layers. The first layer is a stagnant layer of hydrated ions (Helmholtz model) and the second layer is a diffuse layer of the ionic species (GC model). Figure 5 shows the reproduced sketch of the Stern EDL model from his paper [41]. According to his model, the electric potential changes linearly from the solid surface (ψ 0 or ψ s ) to a Stern plane potential (ψ 1 or ψ d ). The ions inside the Stern layer are assumed to be immobile, while the ions beyond the Stern layer are assumed to be point charges (according to the GC-EDL model). It is worth noting that the distribution in the diffuse layer of the Stern model is determined via the Boltzmann distribution equation.
The boundary between the immobile (Stern layer) and mobile (diffuse layer) layer in the Stern model is assumed to be the shear plane in which the zeta potential (ζ ) represents the electric potential on this plane. The real position of the shear plane has been the subject of long debates, which is out of scope of the present work [17,[42][43][44][45][46].
Although the Stern EDL model is considered a significant step forward in our understanding of the surface charge at the solid-liquid interface, it still deficient in several aspects, such as ignoring the chemical reactions that take place on the solid-liquid interface, the solution pH effect, and, importantly, it does not provide a method to calculate the surface charge and zeta potential based on the solution pH and ion concentration.
Approximately 50 years after Stern, Davis et al. [47] attempted to explain and determine the surface charge and electric potential on the Stern EDL model by proposing a detailed surface complexation model for oxide surfaces. They showed that the surface charge is dominated by surface complex formation of ionizable sites and electrolyte ions. However, for a very dilute solution, the surface ionization is significantly related to the concentration of the hydronium (H + ) and hydroxyl (OH − ). The remarkable conclusion of their model is that the electrolyte ions and the hydronium and hydroxyl concentrations are working jointly to determine the surface charge and electric potential via chemical reaction with the surface sites. Their novel results proposed that adsorption density, surface charge, and zeta potential could be estimated simultaneously. The surface of the chemically active materials, such as quartz, clay, etc. not only could be simply ionized but also the surface complexes formed because of the association and dissociation of the hydronium and hydroxyl ions from the solid surface [48].
In this model, it is assumed that the Stern layer consists of two layers, both belonging to the immobile part of the EDL. The complex of the Stern layer with the diffuse layer is called the electrical triple layer (ETL) model. Figure 6 depicts the configuration of the ETL model, in which the Stern layer is divided into two layers via the inner-and outer-Helmholtz planes (OHPs) following the Stern-Grahame electrostatic capacitor model of the inner region [45,48]. In this model, it is assumed that the OHP coincides with the starting edge of the diffuse layer or, in other words, the shear plane.
According to the model that Davis [47] proposed, which was later developed by Kitamura et al. [15], the contribution of the salt-ion adsorption to the surface charge on mineral surfaces is based on the following four chemical reactions [49]: Figure 6. Configuration of the ETL model. In this model, it is assumed that the ions could be distributed in three layers that lie on two planes near the charged surface. The 0-plane, β-plane, and dplane are the inner-Helmholtz plane, outer-Helmholtz plane, and the starting edge of the diffuse layer, respectively [47].
where K int in each equation represents the associated equilibrium constants for the reactions. Equations (10) and (11) explain the impact of the hydronium concentration. If the concentration of H + increases, the chemical reactions (10) and (11) will move to the left side. Similarly, the impact of the metal ionic species (cation M + and anion A − ) is explained by Eqs. (12) and (13). It is noteworthy that in the ETL model, the protonation rate of the siloxane group is very low for the pH range of 3-9. As a result, we can ignore the fourth reaction (Eq. (13)) in the set of ETL equations described above. In the ETL model, the surface charge and electric potential on the three planes are calculated by considering the continuity equation for the surface charge density on the solid surface, the Grahame equation [45], and the differential capacitance of the inner-and outer-Helmholtz layers. Although the ETL model proved to be a great step forward in figuring out the surface charge and electric potential at the chemically active solid surface and the aqueous solution interface based on the properties of the solid surface, it cannot explain some experimental observations such as the impact of solution temperature on the zeta potential [50] or the concentration-dependent electrical conductivity of ultranarrow channels (approximately 2 nm) [17]. With this aim, Alizadeh et al. [17] further developed the ETL model in which an extra layer was added between the zeta potential (ZP) plane and the OHP, called the buffer layer (BL) (Fig. 7). The significant property of this layer is that its position from the solid surface has a flexible nature, which is a function of the bulk ion concentration and solution temperature. In their work, Figure 7. Electrical quad-layer model, in which a new layer was added between the zeta potential plane and the OHP, namely the BL. The figure has been reprinted from [17]. The positions of different planes are shown by X and the thickness of each layer is represented by δ.
they showed that the zeta potential plane (ZP) recedes to the bulk solution or approaches the solid surface when the bulk ion concentration decreases or increases, respectively.
The same chemical reactions were considered as mentioned above for ETL (Eqs. (10)- (13)). The BL properties were integrated into the model by assuming a new capacitance, which was defined as where C BL , Q ZP , and ψ represent BL capacitance, the surface charge on the starting edge of the diffuse layer, and the electric potential on the OHP and ZP, respectively. It is worth noting that the BL capacitance is where ε 0 ε r is the aqueous solution electrical permittivity; and δ BL is the BL thickness, which is a fitting parameter. To determine the thickness of BL, one needs to employ the ab-initio methods (i.e., molecular dynamics [MD] simulations) to gain a deeper understanding of this layer. However, as the aim of Alizadeh et al. [17] work was to present a macroscopic model, the thickness of BL was defined as a fitting parameter. They demonstrated that the electrical quad-layer model could predict the measured zeta potentials and surface charge (Fig. 8) as a function of bulk solution properties (i.e., bulk ionic concentration and solution pH). The intense electric field within the electric double layer could significantly alter the physical properties of the water solvent from its bulk values. Specifically, the orientation of water molecules in the vicinity of the charged surface favors hydrogen bonding interactions, thereby increasing the shear force in the flow direction normal to the surface. From an activation energy perspective, the imposed electric field de-creases the vibration frequency of molecules giving rise to an increase in the activity energy E a and, as such, an increase in the apparent viscosity (according to Reynolds' equation for the viscosity of liquids μ, where the coefficient A is constant): In 1939, Andrade and Dodd were the first to experimentally demonstrate this phenomenon for polar solvents by measuring the flow rate of a pressure-driven flow through a slit channel with an applied electric field normal to the flow direction [51,52]. They observed that in the case of polar solvents, the apparent viscosity increased considerably with the magnitude of the applied electric field whereas it remained constant for nonpolar solvents. In a subsequent study published in 1950, they quantified this phenomenon and revealed that the increase in viscosity is a function of the applied magnitude of electric field, as seen in Eq. (17) [53]: where μ 0 is the viscosity in the absence of an electric field, E is the local electric field, and f = E a /(k B T ) is referred to as the viscoelectric coefficient, in which E a is the activation energy variation due to the presence of E. At low applied electric fields, that is, f |E| 2 1, the higher order terms in the Maclaurin expansion of the exponential function in Eq. (17) become negligible, resulting in a simplified viscoelectric equation, which shows that the increase of viscosity is proportional to the square of the local electric field: According to Debye and Onsager's theories, f is subject to the degree of polarization of the solvent molecules, which can be expressed as where τ is the dipole moment. In 1961, Lyklema and Overbeek analyzed f for water based on a polarization theory of spherical molecules [42]. Using an average value of structural coefficients from previous data of organic polar solvents, f for water was estimated to be 10.2 × 10 −16 m 2 V −2 . This theoretical estimate was later verified by Hunter and Leyendekkers in 1978, who experimentally obtained a slightly lower f = (5∼10) × 10 −16 m 2 V −2 for water between two parallel sheets [54]. Hsu et al. [18] considered viscoelectric effects in EDL models via computational simulations, following the viscoelectic coefficient estimated by Lyklema and Overbeek, and quantitatively obtained experimental measurements of four classical electrokinetic phenomena: electrophoresis, electroosmosis, streaming potential, and streaming current [18]. They reached a consistent conclusion for streaming potential and streaming current: the measured electrokinetic quantities are independent of the surface charge (potential) as the  and ρ e denote the velocity, electric field in EDL, and space charge density, respectively. The figure has been reprinted from [18]. electric field at the interface is higher than a critical value. As seen in Fig. 9, as opposed to models that assume a constant viscosity (e.g., the basic Stern (BS) model, which predicts that the streaming current will increase with the increase in surface charge), the viscoelectric double layer model (VEDL) shows that the streaming current is insensitive to the surface conditions at high surface charge densities due to the presence of a viscoelectric immobile layer. This implies that when the surface charge is high, the surface charge conditions may not be accurately measured by electrokinetic methods.
Using MD simulations, Qiao et al. pointed out that the increase in viscosity is significant at the interface of a charged surface and an electrolyte aqueous solution [55], consistent with the viscoelectric theory. On this basis, Hsu et al. investigated viscoelectric effects in nanochannels using a continuum VEDL model and indicated that viscoelectric effects are dominant over other effects such as ionic steric effects and dielectric saturation effects [56]. Due to the presence of the viscoelectric immobile layer, whose length scale becomes comparable to the channel dimension, the electroosmotic mobility can even decrease as surface charge increases, consistent with Qiao et al.'s MD analysis. Both MD and continuum VEDL simulations show that enhancement of viscosity directly suppresses ion diffusivities D following the Stokes-Einstein equation: where r denotes the ionic radius. In a separate study, Kaji et al. experimentally observed an increase in viscosity and a decrease in DNA diffusivity in nanospace that supports this behavior [57]. Abundant scientific evidence has shown that viscoelectric effects are of crucial importance in various electrokinetic systems, especially at the nanoscale. However, they are largely overlooked in current literature and deserve more attention to elucidate electroosmotic behavior from a fundamental point of view. Thus far, we discussed the various EDL models that have been developed until now. However, we must note that there is another phenomenon, electroviscous effects, which must not be confused with the viscoelectric effect. The electroviscous effect is essentially a result of driving the ionic species from one side through a charged micro/nanochannel to the other side using an applied pressure gradient. In this phenomenon, EDL will trap some counter-ions to drive from the inlet reservoir to the outlet. As a result, we will have nonzero net charge density at the outlet reservoirs, which gives rise to an electric potential difference between the inlet and outlet reservoirs, called streaming potential [58][59][60]. This nonzero potential difference will generate a backward fluid flow to the pressuredriven flow as the movement of the ionic species pull the water molecules along with them due to the friction forces between ionic species and water molecules [58,61,62]. The net effect of these two driving forces will decrease the flow rate in the pressure-driven flow direction. This reduction of flow rate due to the resistance of the streaming potential is equal to the increased viscosity of the aqueous solution. For further detail regarding the impact of electroviscous effect on the pressure-driven flow in micro/nanochannel, refer to Dongqing Li's textbook [58].
Here it is worth pointing out that the distribution of the ionic species at the vicinity of the charged solid surface will be strictly defined by the volume of counter-ionic species when the solution is concentrated. This phenomenon which is the so-called Steric effect [63,64] has been proposed to include the steric replusion. Borukhov and Andelman [63] showed that the Steric effect could be incorporated into the Posisson-Boltzman equation to consider the ionic volume effects which results in ion volume-dependent counter-ion concentration.
It has been shown that for larger sizes of ionic species, the concentration of counter-ion species at the vicinity of the charged solid-liquid interface decreases significantly. The impact of ion size on the electrokinetic flow has been investigated in several previous works [22,[65][66][67]. For instance, considering the transport phenomena through nanoporous membranes, the interplay of ultranarrow pore and ionic volume size will determine the electrokinetic transport phenomena [67].

Mechanism of EOF pumping
In this section, we are going to discuss how the interplay of the applied external electric field and EDL results in pumping of the aqueous solution through the conducting media. Therefore, we focus on a simple flow through two parallel electrically charged plates, called slit microchannel. Figure 10 shows a schematic 2D illustration of the slit microchannel, which is assumed to be negatively charged. EDLs in the vicinity of two walls are shown in gray and are positively charged due to higher counter-ion concentration.
Recalling what we discussed for Eq.
(1), we should point out that the electrical potential at each point in a straight channel without any entrance effect can be considered as a superposition of the applied external electric field and the EDL's electric field. Consequently, we can write the total electric potential as where φ 0 is the anode's potential, x the length from anode, and E x = φ 0 Lc (L c denotes the distance between the two electrodes, which can be considered as the microchannel length). If we introduce Eq. (21) into Poisson's equation, the electric potential due to the EDL will remain while the external electric field is a linear function of x.
By applying an external electric field, the whole aqueous solution feels a body force equal to ρ e ∂φ/∂x. Evidently, the body force can be zero if the concentration of the counter-and co-ions are equal ( ρ e = 0) or nonzero where the local electroneutrality is violated (ρ e = 0). Therefore, if the zeta potential on the microchannel's walls does not change axially, the Lorentz force because of the external electric field will be the only driving force for the whole solution. Considering u EOF in Fig. 10, the shear force at the middle of the microchannel in EOF can be zero ( ∂u ∂y = 0). This phenomenon is owing to zero Lorentz in the middle of the microchannel ( ρ e = 0). Later, we will discuss the distribution of the electric potential and velocity inside micro/nanochannels and porous media.
Here it is worth pointing out that although this tutorial focuses on EOF through media composed of nonconductive materials, EOF can still be applicable when the solution is in contact with a conductive medium where the surface is electrically polarized owing to the external electric field, known as induced-charge electroosmosis (ICEO) [68,69]. Bazant and Ben [70] theoretically showed that the ICEO can be controlled by a 3D AC applied electric field, enhancing the flow rate of EOF upto 20 times compared to planar AC electroosmosis (AC-EO) pumps. Regarding this category of EOF, we direct the readers to another review for details [71].

EOF through microscale pores and channels
In the previous section, we discussed the history and the physics underlying EOF. In this section, we will focus on the broad applications of EOF in the microscale domain, which consists of channels or the porous media. The use of EOF in microfluidic applications is broadening day by day-from lab-on-a-chip devices for biomedical applications to manipulating fluid flows for logical parts of microfluidic chips such as micromixers and microvalves.

EOF through straight microchannels
Following what we discussed in Section 2, if we fabricate a microchannel from materials with functional groups (i.e., ̶ OH) and connect the microchannel to two reservoirs filled with aqueous solutions, upon applying an external electric field through insertion of two electrodes to the reservoirs, the free electric charges in the vicinity of the solid-liquid interfaces (c.f. Section 2.2) start to move toward the counter-charge electrode. The motion of these ions pulls the solvent molecules due to the frictional forces between molecules.
There are various types of species-driven methods in lab-on-a-chip devices, namely, pressure-driven [72], electrokinetic-driven [73], droplet-driven [74,75], and capillary-driven [76]. Among these species-driven methods, we will focus on the electrokinetic-driven method. EOF is one of the categories of species transport methods under the electrokinetic transport phenomena. As we mentioned above, the main advantages of EOF include pumping of solutions without the aid of mechanical parts and easy manipulation (or, in other words, active control) of the fluid flow. In this section, we will try to demonstrate some of these applications and discuss the underlying mechanism.
Considering the structure of lab-on-a-chip devices, which are usually designed to analyze a sample of the analyte [77], one of the key parts of these devices is the micromixer. The rapid mixing of distinct types of species is vital for biochemistry analysis, synthesis of nucleic acids, and even drug delivery. Nguyen and Wu [78] categorized micromixers into passive and active types. The electrokinetic mixing methods are classified as an active mixing method in which the mixing efficiency can be controlled actively through external manipulation. Electrokinetic mixing has attracted the attention of researchers as both pumping and disturbance (which causes the mixing of the species) can be performed simultaneously. It is worth noting that electrokinetic disturbance can be combined with other species-driven methods, such as pressuregradient method, wherein an AC electric field is usually applied perpendicular or alongside the fluid flow. This method will be discussed in detail later.
As one of the pioneering works in this field, Jacobson et al. [79] demonstrated electrokineticallydriven mixing of species in T-shaped microfluidic channels. In their experiments, the samples (which must be mixed with buffer solutions) were driven by applying external electric fields to the reservoirs. The external electric field drove both the sample and buffer from separate channels (Fig. 11).
Applying an external electric field at the reservoirs of the sample and buffer initiated their respective flows to the intersection. Consequently, we can consider the ionic current through the analysis channel as a function of the ionic current through the sample (I s ) and buffer (I b ) channels based on the principle of conservation of ionic current: In this experiment, the microchannels were designed with the same cross-sectional areas to ensure that the impact of the cross-sectional area on the fluid flow of different liquids could be ignored. It is worth pointing out that the applied electric fields on both sample and buffer are the same. Based on Figure 11. Schematic illustration of the T-type micromixer with electrokinetically-driven samples and buffers [79]. Arrows demonstrate the direction of EOF.
these considerations, the ionic current through each channel can be defined as where W , H, ρ e , and u denote the width and height of the channel, the net electric charge density, and the fluid flow velocity, respectively. Thus, we can define the sample fraction as The sample fraction could be also obtained based on the channel length L using It is noteworthy that the above sample fraction relations are valid when the acquired surface charge on the solid-liquid interfaces could be considered constant. However, this is a simplified assumption that is not realistic in many practical applications, especially when there are asymmetrical bulk properties (i.e., concentration gradient, temperature gradient, etc.).
Jacobson et al. [79] proposed microchips for parallel and serial electrokinetic mixing as well ( Fig. 12A and B). In the parallel mixing design, several parallel reservoirs for samples and buffers were designed with connected T-intersections. By applying external electric fields, the fluids from different reservoirs were made to flow and meet at the intersections. To visualize the mixing procedure, a fluorescence solution was added to the sample solutions and the intensity of the fluorescence signal was considered as an index for the mixing capability of the microchip. In parallel mixing, the T-intersections play a key role. In Fig. 12A, the labels "S," "B," and "A" refer to the sample, buffer, and analysis channels. The numerous T-intersections seen in Fig. 12 enhance the mixing of the species. In serial micromixers, the sample and buffer reservoirs were connected to a serial branch of the microchannels, which increases the meeting intersections of the sample and buffers to increase the mixing efficiency of the species ( Fig. 12C and D).
In this method, a four-way intersection plays a key role in enhancing the mixing of the species (Fig. 13). Similar to the parallel mixing method, we can simply derive the conservation of ionic current for the four-way intersection as The above-mentioned EOF micromixer drives the fluid flow by applying a DC electric field and relies on the geometry of the designed microchips. Although this method is categorized as active micromixing [78], we argue that this should be categorized as "semiactive" micromixing method instead. The main reason in favor of this argument is that the mixing relies enormously on the design of the intersections as well as the number of connected microchannels. However, we must mention that the mixing of the species can still be controlled via an applied DC electric field, which controls the speed of flow of the solutions.
Following the idea that there is a direct relation between perturbation and mixing efficiency, Oddy et al. [80] proposed an interesting design that employed an AC electric field instead of the DC electric field for stronger perturbations in a chamber, and therefore, showed better mixing efficiency. The AC electric field induces oscillating EOF in the microchannels. To determine the oscillating EOF, one may consider the Navier-Stokes (NS) equations under the following assumptions: 1D, low-Reynolds number flow (consequently ignoring the advective effects on EDL), which can be simplified as where u, t, and y denote EOF velocity, time, and normal direction to the microchannel's walls. The AC electric field affects the solution located in the vicinity of the solid-liquid interface (i.e., EDL). This effect can be introduced in the NS equations by assuming a slip boundary condition as In this equation, U HS defines the reference velocity for EOF or the Helmholtz-Smoluchowski (HS) velocity [21], which is the result of solving the NS equation for thin EDL (κH 1, see Eq. (5)) with the nonslip boundary conditions on the microchannel's walls ( u y = 0 = 0 and u y = h = 0). The original steady-state 1D NS equation in the presence of an applied external electric field is written as where the last term on the right-hand side of Eq. (29) represents the applied external electric body force. As we mentioned in the previous sections, the net electric charge density (ρ e ) can be obtained by solving the Poisson-Boltzmann equation and recalling the Debye and Hückel assumption: By introducing Eq. (8) into Eq. (30), we obtain An analytical solution to Eq. (32) can be simply derived as by considering the boundary conditions we mentioned above. Next, let us recall the assumption that the EDL is very thin in comparison with the channel height (κH 1). As a result, Eq. (33) can be simplified as and for y = 0, we have u y = 0 = U HS = − ε 0 εr Ex μ ζ . This is called the HS velocity or the slip velocity at the outer edge of EDL. An analytical solution to Eq. (28) is obtained as [80] u (y, t) = f (y) exp (iωt) .
The idea behind the AC-EOF micromixer is simple and effective. A high-voltage amplifier is connected to one reservoir of a T-type microchip made of polydimethylsiloxane (PDMS) and the other end is grounded. The two samples that must be mixed are introduced through the two other reservoirs using a pressure gradient (Fig. 14). The two samples entering the horizontal microchannel will be disturbed due to the applied AC electric field. The oscillation of EOF at the  [80]. The two reservoirs on the right-and left-hand side of the microchip are subjected to an AC electric field while samples A and B are pumped into the vertical channel. Applying an AC electric field will disturb the two samples that are entering the straight channel and will enhance the mixing of the species. horizontal microchannel will enhance the mixing of the samples (Fig. 15).
Thus far, we have shown that the mixing of species could be enhanced by perturbing EOF with the geometrical design or the AC electric field. In addition to these methods, Alizadeh et al. [81] have shown that the mixing of the species Figure 16. Schematic illustration of micromixers with temperature-patterned walls. The red blocks on the microchannel walls represent high-temperature patterns with nonzero surface charge. The other parts of the microchannel were both kept at inlet solution temperature and zero surface charge. The figure has been reprinted from [81]. could be enhanced by applying temperature-patterned walls. In this mixing design, patterns of high-temperature plates were assumed on the microchannel walls in symmetric and asymmetric arrangements (Fig. 16). Patterns with temperatures higher than the inlet solution temperature are shown by red blocks on the microchannel walls. To drive the species at the inlet of the microchannel, different external electric fields were applied along with the same pressure gradient. The modified NS equation for incompressible fluid is where ρ (kg/m 3 ) is the density of the electrolyte, p (Pa) is the fluid pressure, ν (m 2 /s) the kinetic viscosity, and F (N/m 3 ) is the body force density. It is worth noting that the body force could include all applied body forces such as the electrical body force or the pressure gradients. In this study, the electric body force was defined as F = F e + F p = −ρ e (∇ϕ + ∇ψ ) + ∇P, where ∇ϕ is the external electric field and ∇ψ denotes the internal electric field owing to the distribution of the ionic species.
As the temperature field and the distribution of the electric potential on the microchannel walls are nonuniform, one must solve the Nernst-Planck (NP) equation to obtain the distribution of the ionic and sample species. However, the conventional NP equation is based on isothermal solutions. Therefore, we need to modify it for the nonisothermal scenario. Alizadeh et al. [81] modified the NP equation as where the last term on the right-hand side of Eq. (37) is responsible for the temperature gradient effect. This term represents the thermoelectrochemical migration phenomenon, which is a contribution of both temperature gradient and internal electric field. For a detailed discussion regarding this phenomenon and its impact on ionic distribution, see Alizadeh et al. [81].
The two other governing equations that must be solved are the energy equation and the species concentration equations: and where T (K) denotes the solution temperature, k the thermal conductivity, and c p is the specific heat capacity. In Eq. (39), C s represents the species concentration and D s is the species diffusion coefficient. For mixing the species in this microchannel, we need to solve the NS, NP, energy, advection-diffusion, and the Poisson (Eq. (1)) equations in an iterative coupled numerical procedure. The numerical method that Alizadeh et al. [81] used was the Lattice Boltzmann method (LBM), which is a rather novel computational fluid dynamics method that has drawn considerable attention recently [82][83][84] owing to its capability and simplicity in modeling multiphysicochemical transport phenomena through structured and unstructured media with micro/nanoscale characteristic length.
The modeling results demonstrated that by increasing the temperature difference between the solution and the patterns on the microchannel's walls, the induced vortices influence a major part of the microchannel. Consequently, the mixing of the species is enhanced (Fig. 17).
Thus far, we have discussed the electrolyte flow through charged media triggered by an external electric field. However, in some cases the electrolyte flow is generated by applying both pressure gradient and external electric field [85][86][87][88].
The governing equations for the combined EOF/pressuredriven flow are the Poisson equation (Eq. (1)) and the modified NS equation for incompressible flow (Eq. (36)). Equation (36) can be simplified to a 2D system when straight microchannels with the height (h) is much smaller than the channel width (w). Based upon the aforementioned assumptions, we can derive a modified NS equation as If we introduce the Poisson equation (Eq. (1)) into Eq. (40), then we have Equation (41) is a linear partial differential equation which justifies the decomposition of the total velocity into the EOF and pressure-driven velocities as where u p represents the pressure-driven fluid velocity. Equation (41) reduces to the HS velocity (Eq. (34)) when there is no applied pressure gradient (i.e., dp/dx = 0). If we nondimensionalize Eq. (41) withū = u/u HS ,p = p/( μu HS h ),ȳ = y/h, ψ = ψ/ζ , andx = x/h then we have [85] Considering the superposition principle for linear partial differential equations, Eq. (43) can be solved as Figure 18 shows the nondimensionalized combined EOF/pressure-driven flow velocity (ū) versus the crosssection of the microchannel (ȳ) for different applied pressure gradients. Since the EDL thickness is very thin, for zero pressure gradient (i.e., dp/dx = 0), it is demonstrated that the flow velocity decays fast by getting far from the microchannel's walls and a uniform plug-like velocity develops. However, by applying nonzero pressure gradients, the flow velocity departs from plug-like velocity and which depends on the applied pressure gradient direction. Figure 18. The nondimensionalized combined EOF/pressuredriven velocity along the cross-section of the microchannel for different applied pressure gradients. The amounts of dp/dȳ curves' labels. The figure has been reprinted from [85].
Electroosmosis can not only drive Newtonian fluids but also non-Newtonian fluids such as power-law and viscoelastic [89][90][91] fluids. Generally, the non-Newtonian fluids are defined as fluids in which there is a nonlinear relation between the variation of velocity and the shear stress. EOF of viscoelastic fluid between two parallel plates could be described via the Phan-Thien-Tanner (PTT) model [92] which is a model based on the network theory for polymeric fluids. In the viscoelastic EOF, Dhinakaran et al. [91] utilized the Gordon-Schowalter convected derivative which results in nonzero second normal stress difference in pure shear flow. The impact of EDL and the distribution of the ionic species are considered by solving the Poisson-Boltzmann equation. The governing equations for the viscoelastic flow could be considered as the continuity and the modified Cauchy equations as [91] ∇ · u = 0, where η s is the Newtonian solvent viscosity and τ is the polymeric contribution to the extra-stress tensor which the solvent viscosity is assumed to be negligible in comparison with the polymeric contribution (i.e., η s = 0). The external body force is considered to be obtained by F = ρ e E. Considering the PTT model to take into account the viscoelastic behavior of the fluid, one has a relation between the extra-stress tensor and solvent viscosity as [91,92] where D = (∇u T + ∇u)/2 is the rate of deformation tensor, λ the relaxation time, η is the polymer viscosity coefficient, andτ represents the Gordon-Schowalter convected derivative of the stress tensor defined as where ξ denotes the slip between the molecular network and the continuum medium [92]. Dhinakaran et al. [91] proposed an analytical solution for the set of governing equations. Their solution showed that the normal and shear stresses are approximately zero at the centerline of the channel while they rapidly rise by approaching toward the channel's walls.

EOF through microporous media
In the last section, we discussed the solution transport through the electrically charged microchannels with the application in micromixers as an essential part of the lab-on-achip devices. In this section, we will study the applicability of EOF through complex electrically charged domains such as porous media. Porous media can be defined as any naturally (i.e., underground soil or brain) or artificially (i.e., filtration or fuel cell membranes) produced complex material with a very high surface-to-volume ratio. Applications of such porous media are very broad-from biology [93,94] and environmental science [25,26,[95][96][97][98] to geology [99][100][101][102][103]. This section aims to highlight EOF through microporous media by focusing on relevant applications. Here, we must point out that microporous media are characterized as porous media with an average void size on the order of a few micrometers. Historically, Paul et al. [104] reported that the electrokinetic phenomena could induce high pressures that are suitable for pumping liquids through microporous media. The experiments were performed for a pack of micron-size silica beads (Fig. 19). The difference in the work done by Paul et al. [104] with what we mentioned in Section 2 is in the capability of EOF in generating high pressures through porous structures. Paul et al. have shown that EOF could generate pressures as high as 55 MPa. To consider the porous media, they assumed that a porous medium is a composition of microcapillaries. The characteristic properties of any porous media could be defined by several parameters such as porosity (ϕ) and tortuosity (τ ). Porosity is defined as the volume fraction of the voids and tortuosity is simply defined as the length of the tubules per unit length of the medium. A tortuous porous medium will reduce EOF velocity. Consequently, the averaged EOF velocity in a porous medium for a given porosity (ϕ) and tortuosity (τ ) is [104] where d denotes the effective pore diameter for the packed beads (which has a direct relation with the beads' diameter) and f is a function that is responsible for any overlapping effect of EDL. Clearly, for the nonoverlapped regimes, the average EOF is the HS velocity. They proposed that for a very thin EDL, the flow rate must be an independent function of the pore diameter and is proportional to the applied external electric field (E = V/L).
If we assume that porous media is a combination of many capillaries, then the averaged EOF velocity combined Figure 19. A schematic illustration of micron-sized silica beads. The silica beads will acquire a negative surface charge due to the chemical reaction with the solution (c.f. Section 2.1). EOF is generated by applying an external electric field to both ends of the microporous media. The negatively charged silica beads will generate a polarized layer of the solution (EDL) and the external electric field will push the solution in the vicinity of the solid surface.
with the applied pressure gradient could be obtained as [104] To obtain the pressure gradient due to the applied external electric field, we need to assume a zero net flow, which means EOF is opposed to a pressure-driven flow. Therefore, we have Equation (50) indicates that the induced pressure gradient owing to the applied external electric field is solely a function of the effective pore diameter and the electrochemical properties of the porous medium.
Paul et al. conducted simple experimental measurements to obtain the induced pressure due to the electroosmotic pumping effect. In their setup, fused silica capillaries were packed with nonporous silica (NPS). The porosity of their medium was 0.32. The packed silica beads were placed in between two fluid reservoirs, which were filled with a mixture of 80:20 of acetonitrile and water buffered with 4 mM aqueous sodium tetraborate and pure water with 1 mM borate. To make sure that the EDLs in the pores of the silica beads pack are not on the overlapped regime, we can simply calculate the EDL thickness by employing its definition κ = ( 2e 2 c b i ε 0 εr k B T ) , which results in 5 to 6 nm. As we mentioned above, the measured porosity of the pack was 0.32, which suggests that the average pore diameter to the Debye length would be 10. This indicates that EDLs were not generally overlapped. The experimental setup consisted a packed capillary with one end open to the ambient pressure while a platinum electrode was introduced to the reservoir. The other end of the packed capillary was fitted into a plastic HPLC "T" junction which the side leg was fitted with a platinum wire electrode. Clearly, the external electric field will be applied to both ends of the packed capillary through the designated electrodes. The remaining third leg of the "T" was designed to attach to different diagnostics. Paul et al. [104] measured the pressure generated through EOF by measuring the compression length of an air gap that is trapped in a long capillary tube that is sealed on one end. Figure 20 depicts the experimental measurements for different micro-sized silica beads and the generated pressure gradient for different bead diameters. In this figure, ODS represents OctaDecyl Silyl (C 18 H 37 Si(CH 3 ) 2 ), which is a hydrocarbon group with a cation exchanger C 18 H 37 Si. The dashed lines represent the predicted values of P EOF / V by Eq. (50) which are scaled to the mean of the 3 μm ODS coated NPS beads.
As we know, the acquired or applied surface charge at the solid-liquid interface plays a key role in EOF. Therefore, it is of utmost importance to study the influence of surface charge (σ ) or zeta potential (ζ ) on the ionic species' transport through the microporous media. Historically, researchers have attempted to (I) prescribe the surface charge and zeta potential as an independent parameter, namely homogeneous surface charge, or (II) obtain the local surface charge by solving the EDL models (i.e., Section 2.1), namely inhomogeneous surface charge. Clearly, the first method brings simplicity to the complex electrokinetic transport phenomena through complex geometry. However, these results may neglect the possibly significant impact of the local surface charge on the species' transport properties.
For homogeneous surface charge, there is a large body of literature that focuses on the theory [105][106][107][108][109][110] and experiment [111,112] of EOF through porous media. The common assumption among the early works is that the EDL is very thin compared with the characteristic pore size. Therefore, we can introduce a slip velocity (the HS velocity) at the solid-liquid interface, which includes the impact of the EDL on the solution.
Starting from a pioneering work, Coelho et al. [113] briefly reviewed the previous attempts until 1996 and they attempted to study the impact of EDL thickness on EOF through porous media with an intermediate range of double layer thickness. Coelho et al. attempted to, first, propose the general electrokinetic equations that govern the transport of the ionic species and the solution through an infinite spatially periodic porous medium. Subsequently, the proposed equations were linearized by introducing ionic potentials as proposed by O'Brian [112,114]. Finally, they proposed the general numerical simulation results and compared the simple cases with the available analytical solutions.
The general electrokinetic equations for porous media are the same as those proposed for the straight channel (see Section 3.1). The flux of the ionic species is a superposition of the diffusion (the first term on the righthand side), electromigration (second term), and convection (third term). Equation (51) is obtained by applying the continuity equation, Regarding the fluid flow, we can simplify the NS equation because fluid flow through compact porous media inherently has a very low Reynold's number. This low Reynold's number justifies the utilization of the Stokes equations with transient specification as [113] where μ represents the dynamic viscosity, p is the pressure, and F is the Lorentz electric body force (which we have defined in Eq. (36)). Essentially, we need to solve Poisson's equation (Eq. (1)) to obtain the distribution of the electric potential owing to EDL and the ionic species as charged particles.
To solve the governing equations, we need to define proper boundary conditions. The boundary conditions for porous media in which the solid-liquid interfaces were assumed to be impenetrable are the same as those we introduced for the straight channels. Hence, we have where n denotes the normal unit vector to the solid-liquid interface and u s is the fluid flow velocity on the solid-liquid interface. The first boundary condition represents the zero ionic fluxes into the solid surface. It is assumed that electric potential ψ on the interface is the Dirichlet type. Here, we should note that the boundary condition for Poisson equation could be considered as Neumann-type if we set the surface charge on the solid-liquid interface instead of the electric potential. The relation between the surface charge and the electric potential on the solid-liquid interface could be obtained by considering the definition of balancing the surface charge with the free net electric charge density in the solution: If we introduce Poisson equation into Eq. (55), then we have where at infinity or in symmetrical boundary conditions, we have dψ dy = 0. Finally, we obtain the more general form of Eq. (57), which is given as Consequently, for Poisson equation boundary conditions, we can also prescribe the surface charge and obtain the electric potential with Neumann-type boundary condition.
To complete the problem formulation, we need to determine the structure of the porous media. In this regard, Coelho et al. [113] considered several ordered and disordered microstructures. The ordered microstructures could be considered as the periodic media in which the locations of the grains and, consequently, the pore domain could be described easily by simple vectors. It is worth noting that the diameter of the grains must be the same. As Fig. 21A shows, the 3D space of the ordered porous media could be defined in R 3 as where [R] = [I 1 , I 2 , I 3 ] and [n] = [n 1 , n 2 , n 3 ] are the three basic vectors that characterize the unit cell of the ordered porous medium and the trio of integers that belong to Z 3 , respectively. Coelho et al. linearized the governing equations by assuming that the deviation of the system from equilibrium could be very small when the applied external electric field is small compared to the EDL electric field. As a result, the concentration of the ionic species, electric potential, pressure gradient, and the velocity could be written as where the superscript "0" denotes the amount of parameter in the equilibrium state. It is worth noting that the fluid flow velocity in equilibrium state must be equal to zero, which represents the state when there is no applied external body force. Coelho et al. [113] introduced the above definition of parameters into the governing equations and solved them numerically. For further details regarding the linearized dimensionless governing equation, refer to their work. The most interesting results from their numerical simulation is in the limit of the thin EDL when the slip velocity at the edge of EDL (U HS ) is a valid assumption. The velocity in the bulk region of the pore (far from the solid-liquid interface where the electric potential due to the walls would be near zero) is a function of the local electric field. This implies that the fluid flow in the pore space of the porous media could be considered as an inviscid flow (a flow in which the viscosity or the shear stresses are equal to zero) and is governed by Laplace's equation. However, for the thicker EDL, they showed that the whole pore space contributes to the drag force on the porous media wall owing to the fluid viscosity. Therefore, it is reasonable to make an analogy between the pressure-driven flow and an EOF with thicker EDL which the major solution in pore will be driven by the applied electric field because of nonzero net electric charge density. We will discuss in detail the contribution of comparable EDL thickness with the characteristic length of the channel or pore in the following sections.
Later, Gupta et al. [115] attempted to model the electroosmosis through solid porous media by considering a high zeta potential. The high zeta potential will not allow us to linearize the right-hand side of the Poisson-Boltzmann equation (Eq. (3)). As a result, Gupta et al. attempted to solve the governing equations numerically for different types of porous media.
In another study, Wang et al. [110] investigated EOF in anisotropic porous media, that is, porous media with inhomogeneous directionality of grains. The porous medium in the study was created by considering arrays of ellipses which were packed in a microchannel with 1 μm height. It was assumed that the microchannel walls and the solidliquid interface of the ellipse grains with electrolyte solution were charged equally, that is, ζ wall = ζ ellipse = −50 mV. The bulk concentration of the ionic species was selected to be n b = 10 −4 M. In order to initiate EOF through the microporous medium, an external electric field with strength of E = 5 KV/m was applied. As Fig. 22 shows, we can make the porous medium isotropic or anisotropic by changing the orientation angle θ . Wang et al. [110] simply picked a and b in such a way to have a porous medium with pore sizes 0.7 μm.
They solved the NS equation combined with the Poisson equation where the distribution of the ionic species was obtained by employing the Boltzmann distribution [110]. The modeling result demonstrated that the EOF rate is a function of the orientation angle and the size of the ellipse's semimajor axis. It has been shown that by increasing the ellipse's semimajor axis, the EOF rate increases. This means that if we decrease the pore size or increase the porosity of the porous medium in the y direction, the volumetric flow rate of the solution will increase (Fig. 23A). This behavior of the EOF rate as a function of the semimajor axis of a pack of ellipses is expected when the orientation angle is assumed to be θ = 0 because the packs of ellipses will behave as several parallel microchannels. Thus, more electrolytes will flow in the external electric field direction and less will flow normal to the electric field. In contrast, we expect that by increasing the Figure 22. The schematic illustration of the porous medium which is a pack of ellipses in a microchannel. The height of the microchannel is 1 μm. The figure has been reprinted from [110]. orientation angle (θ ), the EOF rate decreases. Figure 23B shows the impact of the orientation angle on the EOF rate. It is shown that by increasing the orientation angle, the EOF rate will decrease.
Prior to this work, Wang et al. [116] studied the EOF pumping effect in microporous media by considering the particle size effects (Fig. 24A), external electric field (Fig. 24B), bulk ion concentration (Fig. 24C), and the prescribed zeta potential on the particles (Fig. 24D). The microchannel had a 1 μm height and the structured porous medium was considered to be a pack of spheres (Fig. 25). The modeling results revealed that by increasing the particle size, the EOF rate increased significantly. In this structure, it was assumed that the particle sizes are changing while the porosity is constant. This situation results in a different particle number when we change the particle size.
As a step forward in understanding EOF through microporous media, Wang and Chen [117] utilized a random generation-growth method for reproducing the structure of the porous media with prescribed porosity and grain sizes [118,119]. It was assumed that the charge of the solid-liquid interface is homogeneous and was prescribed as an input boundary condition for Poisson's equation. The governing equations are the same as those we mentioned previously (Eq. (1), (52), and (53)). They utilized LBM to solve the governing equations (Section 3.1).
The 3D randomly generated microporous media (Fig. 26) were subjected to an external electric field (E) and EOF initiated through the negatively charged micropores. Their modeling results showed that EOF permeability κ e =ū /E (m 2 /sV), which is defined as the averaged velocity over the strength of the applied external electric field, has a nonlinear relationship with porosity.
As Fig. 27 demonstrates, by increasing the porosity of the porous medium for ε < 0.5, the electroosmosis permeability increases slightly. However, for porosities beyond 0.5, the electroosmotic permeability increases dramatically with the porosity. This behavior of porous media was shown by Wang et al. [120] to be due to the numerical instability of the lattice Boltzmann scheme they selected for 3D modeling of complex geometry.
In addition, Wang and Chen [117] investigated the influence of bulk ion concentration and prescribed zeta potential on the EOF permeability (Fig. 28 ).The modeling results are significant, in that they demonstrate the importance of electrochemical properties of the solution and the solid-liquid interface in determining EOF through microstructure porous media.
Regarding the impact of bulk ionic concentration, modeling results of Fig. 28 demonstrate that electroosmotic permeability always increases with bulk ion concentration. However, this increase in electroosmotic permeability slows down for higher bulk ion concentrations.
We have mentioned (Section 2.1) that the zeta potential and the surface charge of chemically active solids are a function of the bulk ion concentration and the solution pH. Thus, Figure 24. EOF rate versus (A) particle size, (B) applied external electric field strength, (C) bulk ionic concentration, and (D) the particle zeta potential. The figures has been reprinted from [116]. Figure 25. Microporous medium, which is a pack of spherical particles with a structured distribution. The walls of the microchannel were charged as ζ w = −50 mV while the zeta potential of the particles could be changed as a parameter to study the EOF rate. The external electric field (E) and the generated pressure gradient ( P) are shown in this schematic illustration. The figure has been reprinted from [116].
it is of utmost importance to determine the solid-liquid interface electric potential based on the local solution properties. Following Wang and Chen's [117] study, Wang et al. [120] attempted to study multiphysicochemical transport through randomly generated porous media by considering the zeta potential as a result of the local solution properties. For the electric boundary conditions on the solid-liquid interface, they utilized the BS model, which is based on the dissociation of the silanol groups: Behrens and Grier [121] proposed that zeta potential could be expressed as a function of surface charge density (σ ) and defined as Figure 26. Generated microporous structure using the random generated-growth method for a 60 × 60 × 60 grid system. (A) represents the porous medium structure with porosity 0.3 and (B) porosity 0.6. The structures have been reprinted from [117].  where = 8, pK = 7.5, and C = 2.9 denote the surface density of chargeable sites on the silica surface (nm −2 ), the dissociation equilibrium constant, and the Stern layer phenomenological capacity (F/m 2 ). Clearly, the two Eqs. (62) and (63) are a set of nonlinear equations that must be solved numerically, combined with the Poisson-Boltzmann and NS equations. Wang et al. [120] solved the governing equations by using the 3D LBM method. Their modeling results demonstrate the distribution of the electric potential that is obtained by solving Eq. (62) (Fig. 29A) and EOF velocity vectors (Fig. 29B).
This work could be distinguished from the previous one (Wang and Chen [117]) as it highlighted the fact that the zeta potential cannot be changed as an independent input parameter according to bulk ion concentration and solution pH. Consequently, if we need to study the impact of zeta potential on electroosmotic permeability, we need to change the pH or the bulk ion concentration. In their work, Wang et al. [120] showed that EOF permeability sharply changes as a function of the porosity even for ε < 0.3 (Fig. 30). As we mentioned previously, these results are in clear contradiction with their previous numerical prediction (Wang and Chen [117], Fig. 28). These contradictory results were interpreted to be due to the numerical instability of the former study.
Wang et al. [120] studied the effects of zeta potential (Fig. 31A) or, consequently, solution pH (Fig. 31B) on the electroosmotic permeability. Interestingly, the modeling results indicated that there is no electroosmotic permeability saturation for high zeta potentials. This finding was similar to that of their previous work (Wang and Chen [117]) (Fig. 29B). They showed that even for zeta potentials higher than 100 mV, the logarithm of electroosmotic permeability could be considered to be increasing approximately linearly. They concluded that these results, which are in contradiction with their previous study, could be due to numerical instability and not a physical effect. This increase in zeta potential for a fixed bulk ion concentration is due to the increase in the solution pH. In Fig. 31B, the electroosmotic permeability is demonstrated as function of pH.
As a combination of EOF through a straight microchannel and a porous medium, EOF through a microchannel with   [120]. random roughness appears to be a practical candidate [122]. Wang and Kang [122] generated 3D randomly generated microporous media with random roughness on a microchannel. They set three parameters: number density, total volume fraction, and anisotropy of the roughness elements. Figure 32 shows the microchannel with different roughness. For details regarding roughness generation, refer to Wang and Kang's work [122]. The EOF rate in this work was defined as where u z represents EOF velocity in the applied external electric field direction and A is the area of the microchannel cross-section. Figure 32A shows the normalized EOF rate versus the roughness number density of the microchannel. The roughness number density is defined as where N wall and A wall are the total cell number and the area of the microchannel wall, respectively. Increasing the roughness number for a fixed microchannel walls' area and the total cell numbers was carried out by increasing the roughness distribution probability. As we have shown in Fig. 32D-F, an increase in s d would make the roughness microstructure coarser with less connections. The impact of the roughness number density on the EOF rate is shown in Fig. 33A. It is interesting to note that the EOF rate increases linearly with the logarithm of the number density of roughness. These results indicate that a more connected or dense roughness of the microchannel will decrease EOF pumping effect. In Wang and Kang's work, the roughness volume fraction was defined as the ratio of the total volume of the roughness to the fluid volume for the microchannel with smooth walls. If we fix the number density of the roughness as n R = 360/μm 2 for the squares and n R = 36/μm 2 for circles in Fig. 33B, then the numerical modeling results revealed that the EOF rate generally decreases with increase in the volume fraction for both the number densities of the roughness scenario. To understand the impact of volume fraction on the EOF rate, we can imagine that increasing the volume fraction means that the total volume of the roughness increases for a fixed microchannel height. These results are reasonable as by increasing the volume fraction of the roughness, the effective area for electrokinetic transport is decreased. Moreover, it has been shown that by increasing the number density of roughness, the EOF rate will increase. These results are consistent with the impact of the number density on the EOF rate (Fig. 33A).
The works discussed above (except Wang et al. [120] in which the solution pH determines the zeta potential and surface charge of solid-liquid interface) investigated charged microporous media where the surface charge was prescribed and homogenous. This means that the surface charge was an input parameter and did not change with the solution properties. However, Zhang and Wang [123] investigated electroosmosis in an inhomogeneously charged microporous medium. An inhomogeneously charged porous medium refers to the fact that the surface charge on the solidliquid interface must be obtained based on the local solution properties (i.e., pH, temperature, and bulk ion concentration). Therefore, one of the EDL theories that we introduced in Section 2.1 must be employed. Zhang and Wang [123] Figure 31. Impact of (A) zeta potential when the bulk ion concentration is n b = 1 × 10 −5 M and T = 293 K and (B) solution pH on electroosmotic permeability when ε = 0.14 and n b = 1 × 10 −5 M. The figures have been reprinted from [120].  obtained the local surface charge by using a 1-pK model, which is a simpler solution compared to complicated models such as the triple layer model [15,47,124] or the quad layer model [13,17]. However, the details of the model that was used by Zhang and Wang is not the focus of this tutorial. By following the random generation-growth method, they reconstructed a 3D microporous media as shown in Fig. 34. The external electric field was applied in the x direction and EOF was in this direction. The same governing equations as the previous works were solved, except the boundary condition for the Poisson's equation in which the zeta potential was ψ S = ζ (x, y, z).
As we discussed above as well, one challenging fact about EOF through straight channels or the porous media is to Figure 34. 3D microporous medium that was generated by the random generation-growth method. The black parts represent the solid and the blue parts are the voids filled by the solution. EOF will be in the x direction. This figure has been reprinted from [123].
figure out the overlapping of EDLs. Therefore, it is essential to realize when the overlapping happens. Therefore, Zhang and Wang [123] proposed a parameter inspired by the Knudsen (Kn) number. Generally, the Kn number is defined as the ratio of the mean free path of the molecules to the characteristic length of the space, which could be the channel height or the pore diameter. This number determines whether continuum mechanics or statistical mechanics must be utilized to model the fluid mechanics. Inspired by the Kn number, they proposed the M number, which is equivalent to the ratio of the EDL thickness to the characteristic length of the domain M = λ/L. Based on this number, EOF can be divided into four regimes as ( As Fig. 35 demonstrates, by increasing the EDL thickness, most of the channels' bulk will be charged and influenced by the charge of the solid-liquid interface. It is worth pointing out that any further increase of the EDL thickness, which can be done by decreasing the characteristic length of the domain or decreasing the bulk ionic concentration, will lead to monotonic distribution of the electric potential across the channel or pore at M 1. Therefore, we can assume that the electrical charge of the whole channel or pore is identical to the charge of the solid-liquid interface. To investigate the effect of the applied external electric field and the inhomogeneity of the surface charge on EOF, the averaged EOF velocity for different applied electric field strength and two pH gradients were studied for a porous media with porosity ε = 0.46 (Fig. 36). The modeling results show that for E < 20 V/m, EOF shows a nonlinear behavior with respect to the applied electric field. This effect could be interpreted by considering the fact that the diffusivity of the hydrogen (H + ) and the hydroxyl (OH − ) ions are different. When we increase the applied external electric field, the movement of hydrogen ions toward the outlet will be greater, and the local distribution of the pH will be affected by the applied external electric field. While the local pH distribution close to the outlet decreases, as we know, the local absolute zeta potential will decrease. The interplay of the applied external electric field and lower influence from the solid-liquid interface at the outlet of the porous media will generate a nonlinear behavior of EOF velocity. However, by increasing the applied external electric field, the impact of the electromigration phenomenon will be dominant. As a result, EOF velocity can be determined by a linear relationship with E (see Fig. 36 for E > 20 V/m).

EOF through nanoscale pores and channels
In the previous section, we discussed EOF and its application when the characteristic length of the domain is on the order of a few micrometers. This scale explained how the solution could be transferred through nonoverlapped EDL regimes. We introduced different practical applications, including micromixing, EOF through porous media with applications in energy (e.g., enhanced oil recovery) and environment (underground water remediation), as well as two methods to characterize the electric charge at the solid-liquid interface. The first one, which is generally employed for the sake of simplicity, prescribes the surface charge as a spatially nonvariable parameter while the second one obtains the local surface charge based on the local solution properties. Clearly, the second strategy is more realistic in considering the effects of local solution property on EOF in the presence of applied pH or concentration gradients. In this section, we introduce EOF through straight channels and porous media for characteristic lengths on the order of a few nanometers. At this scale, the transport medium guarantees overlapping of the EDLs, due to which interesting phenomena begin to emerge.

EOF through a straight nanochannel
As briefly alluded to earlier, we are interested in studying EOF through nanometer sized channels due to the strong overlapping of the EDLs. This results in a nonuniform velocity distribution along the nanochannel's height due to the strong interaction of the major solution with the applied external electric field [125]. There is a large body of literature that investigates the physics underlying the transport of ionic species through confined domains such as nanofluidic channels [125][126][127][128][129][130][131][132][133][134]. However, we aim to focus on the fundamentals of electroosmosis and the impact of the solid-liquid surface charge on EOF.
Historically, Burgreen and Nakache [2], for the first time, developed the theory of electrokinetic flow in ultrafine capillary slits. They extended the general theory of the electrokinetic flow, which was earlier limited to the nonoverlapped regime of EDLs. Later, Levine et al. [135] extended Burgreen and Nakache's work in which the charge regulation phenomenon due to the overlapping of EDLs was taken into Figure 35. EOFs based on the EDL overlapping regimes. The blue lines represent the distribution of the electric potential owing to the charged solid-liquid interface. The figure has been reprinted from [123]. Figure 36. EOF velocity through the porous media with porosity 0.46. To induce surface charge inhomogeneity, two pH gradients were employed in which the black symbol line demonstrates the pH from inlet to outlet identical to 6 and 8 and the red symbol-line demonstrates the pH from inlet to outlet equal to 5 and 9. These results have been reprinted from [123].
account. Generally, charge regulation is attributed to the surface charge variation owing to reduction in distance between two charged plates [136]. This phenomenon is usually observed in overlapped EDL regimes where the conduction current exceeds the conduction in the bulk solution [137]. Furthermore, Levine et al. [135] considered the convection current generated due to the applied external electric field. Following Levine et al. method, we are going to introduce the semianalytical solution for EOF through a confined channel, which is considered to be the space between two parallel symmetrically charged infinitely big flat plates (Fig. 37).
EOF will be investigated as a function of the applied external electric field in the x direction, which is parallel to the flat plate surface. It is assumed that the solution is under noslip boundary conditions on the surface of the plates (y = 0). Considering steady-state conditions, Eq. (53) reduces to where due to the symmetry of the system in the y direction, Eq. (67) could be solved for the half-domain as 0 ≤ y ≤ h. The boundary conditions for the NS equation (Eq. (67)) are Here similar to what we did for the microchannel, the impact of external electric field and the EDL electric field is decoupled, which gives rise to where L p denotes the length of the plate in the x direction. It must be noted that this decoupling of the electric potential is only valid when the applied external electric field is small [135] or when ignoring the entrance effects. This assumption will let us consider the total electric potential as a superposition of the applied external electric field and the electric potential owing to the electrically charged solid-liquid interface. Consequently, we can easily solve Poisson's equation for ψ as well. The boundary conditions for Poisson's equation are proposed as If we substitute Eq. (1) into Eq. (67), then we have If we introduce the velocity and the electric potential boundary conditions into Eq. (71), then we have an analytical solution for EOF velocity as where P = − dp dx denotes the applied pressure gradient. Equation (72) indicates that the total electrolyte velocity is a result of both pressure-driven and electroosmotic fluid flow. While we are interested in EOF velocity, we can set P = 0, which Eq. (72) gives rise to It is interesting to note that Eq. (73) reduces to the HS velocity if we assume that 2h λ. Based on this assumption, we can assume that a major part of the channel will not be electrically charged and, as a result, we have ψ (y) ≈ 0. Hence, Eq. (73) reduces into HS velocity U HS = − Eεr ε 0 μ ζ . By introducing U HS into Eq. (73), we have Figure 37. Schematic illustration of the nanochannel, considered to be the space between two parallel flat plates that are equally charged. The external electric field is applied to the x direction. The height of the nanoslit is 2h.
If we are interested in the mean EOF along the nanoslit cross-sections, we can simply integrate both sides of Eq. (74), which gives rise to [2,135] where G is defined as where the function G represents the ratio of the mean electric potential to the electric potential (zeta potential) of the solidliquid interface. Thus far, we introduced an analytical solution for EOF velocity through a confined nanoslit. However, the only remaining challenge is to obtain the distribution of the electric potential along the cross-section of the nanoslit. The Poisson-Boltzmann equation could be utilized for this specific problem, where ψ could be obtained as where ψ * and y * are the dimensionless formats of the electric potential, defined as ψ * = ψ/V T and the unit of length as y * = y/λ, respectively. Equation (77) could be rewritten in the form of Levine and Suddaby [138] demonstrated that there is an analytical solution to the rapidly converging series for the electric potential as where k is defined as k = exp(−ψ * (h)) and K (k) = F ( π 2 , k). F is the incomplete elliptic integral of the first kind and defined as [135] If we introduce Eq. (79) into Eq. (76), then we have a relation for G as where Thus far, we introduced an analytical solution for EOF velocity and the electric potential when the electrolyte is moving through a nanoslit under application of an external electric field. In another attempt to characterize electroosmosis through nanochannels, Pennathur and Santiago [125] proposed an analytical solution for EOF and experimentally validated their theoretical study [126]. Regarding the theoretical work, the NS equations that demonstrate the balance of the viscous flow and the applied Lorentz electric force on the solution were considered to be where μ is the dynamic viscosity and denotes the total electric potential, which could be decoupled (following Probstein [139]) to the external electrical potential due to the applied electric field and the internal electric field because of the EDL effect as = φ + ψ. Integration of Eq. (83) gives rise to Eq. (74). By double integration of Eq. (74) across the cross-section of the nanochannel, Pennathur and Santiago [125] obtained the area-averaged velocity by numerical simulation. However, an analytical solution for low zeta potentials is available based on the Debye-Hückel theory [139]: where Here it is worth reminding that EOF that is proposed by Eq. (84) is only valid for very small zeta potentials when ζ V T .
One way to analyze EOF through the nanofluidic channels is to consider the electrolyte solution with some analyte ion. For instance, let us assume an electrolyte solution with fully ionized background electrolyte with ionic species A and B. The general NP equation with electromigration-diffusionadvection contributions can predict the distribution of all the ionic species. As a low concentration of the analyte ion S does not affect the distribution of background ionic species, we can separate the governing equation for background ionic species A and B from the analyte ion S. With this aim, Pennathur and Santiago [125] developed a relation to explain the analyte ion velocity as where ψ c represents the centerline of the nanochannel electric potential. z S , F, and ν S represent the ionic valence of the sample analyte, Faraday's constant, and the ion mobility, respectively. Equation (86) denotes the average velocity of the sample analytes as a result of EOF flow due to the background solution movement together with the electrophoresis movement of the sample analytes due to the applied external electric field. Figure 38 shows the normalized mobility of the ionic species in a nanochannel with respect to the mobility in microchannel as a function of 2λ/h. It is interesting to note that mobility in both nanochannels with 40 nm and 100 nm height, increases with increase in the thickness of EDL. However, the increase in ionic mobility has a peak at 2λ/h ≈ 0.5 and any further increase in EDL thickness (i.e., stronger overlapping of EDLs) decreases mobility. One reason behind this is that upon an increase in overlapping of EDLs, the first term on the right-hand side of Eq. (86) tends to zero because the electric potential in the bulk solution tends closer to the zeta potential (1 − ψ (y) ζ ) → 0. Consequently, the electrophoretic phenomenon plays a key role in driving the species from the inlet to the outlet of the nanochannel.
The electroosmosis inside the nanochannel could be employed to separate the species with a different electrical charge. As Eq. (86) indicates, the positive sample analytes experience lower velocity in comparison with the negative analytes. This simple fact was used by Garcia et al. [140]   to demonstrate molecular separation in nanoscale fluidic channels. They analyzed the electrokinetic transport of the negatively charged species by introducing a negatively charged dye, Alexa 488, and EOF with a neutral dye, rhodamine B. Their measurements revealed that for a negatively charged nanochannel, the electroosmotic velocity toward the cathode (negative electrode) of the negatively charged dye is higher than the neutral dye. This behavior of the negatively charged species inside a nanochannel is anomalous. It is the opposite of what was observed with microchannels, where the neutral dye was transported faster than the negatively charged dye. The main reason for this phenomenon is that for the microchannel, the electrophoretic drag force will decrease the speed of the negatively charged species while the neutral species move faster. This interesting experiment demonstrated the impact of the transport phenomena in nanochannels with overlapped EDLs and microchannels with nonoverlapped EDL regimes. To explain their experimental measurements, the authors followed the same theory as we discussed above for EOF through nanochannels.
The same relation as Eq. (74) for EOF velocity was proposed while the local electric potential was obtained analytically as [39] ψ * (y) = 4 z tanh −1 tanh ζ * 4 exp (−κy) where ψ * is an approximate solution for the electric potential under three assumptions: (I) 2D flow in (II) channels with parallel walls and (III) weak EDL interaction. By introducing the electric potential (Eq. (87)) into Eq. (74), Garcia et al. [140] obtained the concentration-weighted velocity of a dye as a function of EOF and electrophoretic velocities: whereC = 1 h h ∫ 0 C(y)dy and C (y) = n b exp(−zψ * (y)).
As Fig. 39 demonstrates, on one hand, for the neutral species, there is no peak velocity in the middle of the nanochannel. This phenomenon is expected because the neu-tral species do not feel any electrostatic force from the charged solid-liquid interface. Consequently, there would be a uniform distribution of the neutral dye at the middle of the nanochannel and uniform concentration-weighted velocity at the major part of the nanochannel (Fig. 39, curves (3) and (4)). In contrast, the negatively charged dye will mostly repel the middle of the nanochannel owing to the negative surface charge on the solid-liquid interface, which makes theC/C(y) smaller and consequently increases concentration-weighted velocity at the centerline of the nanochannel (Fig. 39, curves (1) and (2)). Furthermore, Fig. 39 demonstrates the impact of the nanochannel height on the concentration-weighted velocity. It shows that the maximum electrokinetic velocity differences will be in the approximately 60 nm height or when κh ∼ 4. For heights below this number, the interaction of the EDLs will increase and, as a result, an electric potential along the cross-section of the nanochannel will be developed to the zeta potential. This will lead to the development of EOF velocity profiles.
The employment of a nonuniform EOF can enhance the efficiency of protein separation in nanochannels as well as the efficiency of molecule detection in nanopores. Considering a nonuniform zeta potential distribution along the surface of cylindrical capillaries by partially coating with a polymer layer, Herr et al. investigated the induced pressure gradient and nonuniform EOF profiles along the axial direction [141]. Similarly, distribution of surface charge and EDL thickness can be generated by introducing salt concentration or pH gradients in silica channels [142][143][144][145], enabling For the ordered porous media in the vicinity of the perm-selective membrane, the strong EOF (red arrows) and back pressure-driven (green arrows) electrolyte solution will initiate a salt (blue area) depletion region (white area) and the vortices are restricted to the space between grains. (C) For the porous media with a random distribution of the grains, the vortices are generated around the grains. Reprinted with permission from [160]. Copyright (2013) American Chemical Society. Figure 43. Three regimes of ionic species transport through a microchannel, which is dead ended via a perm-selective nanoporous membrane. For further details see [161]. high-performance label-free separation of proteins in nanochannels [146,147]. In terms of nanopore sensing, it has been experimentally demonstrated that both capture rate and translocation time of DNA molecules in a nanopore can be simultaneously enhanced by adding a salt concentration gradient [148]. Hsu and Daiguji theoretically investigated electroosmotic behavior in a cylindrical nanopore channel when a salt concentration gradient exists in the axial direction [149]. Due to the variation in EDL thickness, a nonuniform EOF was induced along the surface, which yielded an induced pressure gradient near the nanopore centerline to satisfy mass conservation, as shown in Fig. 40. As a result, EOF became weaker near the molecule entrance (downstream), facilitating the capture, while the amplified EOF near the trans reservoir in the nanopore (upstream) extended the translocation time. Moreover, it was found that a conical geometry can enlarge this effect, giving rise to an induced reverse EOF (IREOF). This unique characteristic enables flexible control of flow behavior in a nanoconfined space for molecule manipulation in versatile bionanosensing systems by altering the bulk solution conditions. Not only can the concentration difference be directly generated by manipulating the bulk conditions, but a local concentration can be induced at the junctions between the nanochannels/nanopores and the microchannels (or reservoirs) when an electric field is applied over ion-selective channels, due to ion concentration polarization [150]. As seen in Fig. 41A, a clear concentration difference occurs that increases with the magnitude of the applied electric field. The large electric field and concentration gradient across an ultrathin nanopore yield significant ion separation, known as transport-inducedcharge (TIC) [151]. Note that both induced charge and surface charge are negative and, hence, a reversal of EOF appears when increasing the applied electric field, as shown in Fig. 41B. At low applied electric potentials, the EOF behavior is dominated by the charge at the interface (see Fig. 41C) whereas EOF direction along the axis is governed by TIC when the applied electric potential is high (see Fig. 41D).
This nonlinear electroosmotic behavior has an enormous influence on DNA molecule translocation behavior through 2D nanopores, where a threshold voltage is needed for a translocation event to occur to overcome the opposite EOF [152]. An ion concentration difference between the nanopore Figure 45. (A) Measured zeta potential as function of Tris molarity for the packed silica nanospheres. EOF mobility as a function of (B) Tris molarity and (C) the ratio of the effective pore size to the EDL thickness. The experimental measurements carried out for a different pack of solid silica nanospheres with distinct average pore sizes. The results have been reprinted from [162].
inlet and outlet could be induced when an electric field exists within a nanochannel. The phenomenon is known as ion concentration polarization [153]. When a current is established in a negatively charged ion-selective nanopore, cations at the anode end are transported through the nanopore to the cathode end. Anions at the cathode end are repelled from the electrode and accumulate at the nanopore junction. Conversely, ion concentration at the anode end near the nanopore junction decreases. As a result, a solute concentration difference between the nanochannel junctions is established.

EOF through nanoporous media
In the previous section, we introduced the theory of EOF through nanofluidic channels and discussed some relevant applications. We showed that the height of the channel plays a key role in multiphysics transport phenomena and that EOF velocity depends on the ratio of the nanochannel height to the EDL thickness. In this section, we discuss transport phenomena through unstructured media, such as nanoporous membranes and rocks, which has practical applications in water treatment [154][155][156][157], underground remediation of toxic ionic species [26,158], and biological applications [159].
The transport of ionic species through nanoporous media is of great interest theoretically and experimentally. For instance, Deng et al. [160] investigated the flow through a silica glass frit, which is a negatively charged porous medium. In their experimental study, the ionic current was driven from an anode in a reservoir to a perm-selective (or cationselective) membrane, namely Nafion (Fig. 42A). The experimental setup was a sandwich of a reservoir, a 1 mm silica glass frit with an average pore size of 500 (nm), and a cationselective membrane. They stated that the transport of ionic species through an electrically charged nanoporous media is under the impact of the surface charge and EOF ( Fig. 42B and C). Considering surface conductance (SC), Deng et al. [160] showed that the transport of the ionic species is mainly carried out by the surface charge for ultra-confined porous membranes or low salinity electrolytes. This means that for the strong EDL overlapping regime, the surface charge determines the transport of ionic species through the porous media because the pores are electrically charged and makes it a counter-ion-selective medium.
However, for the weak EDL interactions that could occur at larger pore sizes or more concentrated solutions, the transport of the ionic species is mainly carried out via EOF while the conductivity of the bulk solution is considerably higher than the EDL conductivity. This principle was first predicted theoretically by Dydek et al. [161]. They showed that the transport of the ionic species is mainly a function of the connected microchannel to the impermeable membrane. They found that the height (depth) of the connecting microchannel to the impermeable porous membrane specifies the ionic transport regime, which could be the SC (Fig. 43A), EOF (Fig. 43B), or electroosmosis instability (EOI) (Fig. 43C). For further useful information, refer to Dydek et al. [161]. We mentioned this phenomenon to emphasize the interesting phenomena observed at the interface of the nanoporous membrane and a microchannel.
EOF through the nanoporous media could be characterized by introducing a fluorescent buffer solution to examine EOF velocity through nanoporous media with pores on the order of the EDL. With this aim, Bell et al. [162] carried out experimental measurements on EOF velocity for a nanoporous membrane made as a pack of solid silica nanosphere with effective pore sizes from 104 nm down to 8 nm (Fig. 44).
They utilized theoretical EOF mobility, which was defined by Ref. [163] as where μ EOF , φ, ε, μ, ψ s , and χ denote EOF mobility, the volume fraction, electrical permittivity of the solution, dynamic  viscosity, electric potential on the solid surface (which is considered to be the zeta potential), and the correction factor (which is a function of r e f f λ ), respectively. According to Levin et al. [5], for thin EDLs or larger effective pore sizes, χ = 1; this means that EOF mobility will be solely a function of zeta potential. However, for thicker EDLs or smaller effective pore sizes, it can be determined in a way to consider the influence of the EDL overlapping on EOF mobility. Here, we should point out that EOF mobility is defined as the ratio of the average EOF velocity through the nanoporous media to the applied external electric field μ EOF = u EOF E . Levin et al.'s experimental measurements revealed that EOF mobility is not only increased by increasing the Tris concentration but also by the effective radius (Fig. 45B). Moreover, they showed that EOF mobility for all nanoporous media with distinct effective pore sizes will be reduced to a single curve if we plot EOF mobility as a function of r e f f /λ. This behavior of silica nanoporous media means that EOF mobility not only depends on the surface potential but also the overlapping of the EDLs [162]. EOF mobility, which was predicted by Eq. (89), is supported by experimental data Fig. 45C) where zeta potential of the packed silica nanospheres was measured for different packs (Fig. 45A).
Similar to what we discussed for ionic transport through microporous media (Section 3.2), transport through nanoporous media can be investigated by considering the local solution properties and, consequently, the local surface charge impact on EOF. As we discussed above, the experimental study has shown that transport in porous media with an overlapped regime is considerably different from that in Figure 48. Schematic illustration of in situ EKR of groundwater from biological contaminants. The two transport mechanisms shown here are (I) electromigration and (II) electroosmosis. The figure has been reprinted from [165].
nonoverlapped EDL regime that is typically found in microporous media. With this aim, Alizadeh et al. [96] conducted a theoretical study to investigate effect of local solution properties on the macroscopic transport phenomena in tight porous media with an average pore size on the order of a few nanometers. Three-dimensional unstructured porous media (Fig. 46A) were generated via random generation-growth method that we discussed in Section 3.2 with four porosities ε = 0.3, 0.4, 0.5, and 0.6 ( Fig. 46B-E). Multicomponent ionic species were introduced at the inlet of the nanoporous media, which could be transferred from the inlet to outlet by applying an external electric field. To obtain EOF, it is essential to solving the Poisson, NP, and NS equations together. The impact of the local surface charge was incorporated into Poisson's equation as a boundary condition by employing the ETL model (Section 2.1).
Their modeling results showed a non-negligible impact of the local surface charge on different aspects of transport phenomena through the nanoporous media. As the focus of this tutorial is on EOF, we only present EOF velocity along with the porous media (Fig. 47). As shown in Fig. 47, the inhomogeneous charge distribution and the porosity of the porous media significantly affect EOF velocity. It was shown that by increasing the porosity of the nanoporous media, EOF velocity increases for both homogeneous and inhomogeneous surface charge scenarios. However, EOF velocity is higher for the inhomogeneous scenario compared to the homogeneous scenario.
As we mentioned at the beginning of this section, one of the interesting applications of EOF is in driving underground water from one site to another desired site. Underground solute transport via EOF is one of the great advantages of this transport phenomenon. The principle of underground EOF is simple and straightforward. By applying an electric field, which only needs injected electrodes into the desired sites, we can collect the hazardous contaminants (i.e., arsenic) to a specified place and then remove them completely. This method, called EKR, has drawn considerable attention in recent years because of its ability in removing both organic and inorganic contaminants in low-permeability domains such as underground soil [164]. There is a large body of literature that focuses on EKR. Its mechanism is the same as what we discussed for channels and porous media, which are based on (I) electromigration (displacement of charged species) and (II) electroosmosis (movement of a fluid flow owing to the nonzero net movement of the ionic species) (Fig. 48). A mechanistic understanding of the multidimensional transport of organic and inorganic contaminants coupled with the reactions in porous media has drawn attention recently [165][166][167]. In this regard, Sprocati et al. [166] modeled the electrokinetic transport with biogeochemical reactions in porous media for continuum scale by combining COMSOL Multiphysics and PhreeqcRM. The former was used for solving the Poisson-Nernst-Planck equations and the latter [168] for solving the geochemical reactions.

Concluding remarks
EOF has proven to be a practical fluid pumping method for transferring water and charged species through tiny slit channels or pores that our understanding of this unique phenomenon has been developing for more than two centuries. The key role in EOF plays by a charged layer of solution called the EDL that forms in the vicinity of a charged solid surface and consists of ordered water molecules and counter-ions. In addition to the role of EDL, the intermediate domain has also demonstrated an interesting and unique impact on EOF. The relative size of the channel/porous media with respect to the thickness of EDL is also a critical parameter to determine EOF. This is the reason why the transport phenomena of water and ionic species at microscale and nanoscale are largely different. In this tutorial, we summarized the works that consider this effect for both straight channels and porous media.
The rising interest in employing EOF for emerging applications has encouraged scientists to conduct further theoretical and experimental studies. Among the various emerging research directions in this field, an interesting one involves multicomponent multiphase electrokinetic transport with an emphasis on transport in charged micro/nanochannel and porous media with broad applications in energy, biology, and environmental issues. Another interesting direction could be to study the impact of thermodynamical forces because of the concentration, viscosity, and temperature gradient on EOF through micro/nanochannels and porous media by taking into account the impact of these gradients on the local surface charge.

Data availability statement
Data sharing is not applicable to this article as no new data were created or analyzed in this study.