Electrokinetic Control of Viscous Fingering

TL;DR

This paper develops a theory for controlling viscous fingering in electrolytes using electrokinetic effects, predicting that electric current can stabilize or destabilize interfaces in porous media, with applications in oil recovery and soil remediation.

Contribution

It introduces a novel theoretical framework for electrokinetic control of viscous fingering, extending classical fluid instability analysis to include electric fields and electrokinetic phenomena.

Findings

Viscous fingering can be modulated by the electric current to flow rate ratio.

Beyond a critical injection ratio, stability depends only on flow and current directions.

Electrokinetic effects can be used to control fluid interfaces in porous media.

Abstract

We present a theory of the interfacial stability of two immiscible electrolytes under the coupled action of pressure gradients and electric fields in a Hele-Shaw cell or porous medium. Mathematically, our theory describes a phenomenon of "Vector Laplacian Growth", in which the interface moves in response to the gradient of a vector-valued potential function through a generalized mobility tensor. Physically, we extend classical Saffman-Taylor problem to electrolytes by incorporating electrokinetic phenomena. A surprising prediction is that viscous fingering can be controlled by varying the injection ratio of electric current to flow rate. Beyond a critical injection ratio, stability depends only upon the relative direction of flow and current, regardless of the viscosity ratio. Possible applications include porous materials processing, electrically enhanced oil recovery, and electrokinetic remediation of contaminated soils.

Equations

$${\mathbf{u}}^{\pm }=-\frac{b^2}{12{\mu }^{\pm }}\nabla p^{\pm },\nabla \cdot {\mathbf{u}}^{\pm }=0$$

$$[[p]]=\gamma \kappa ,[[\hat{\mathbf{n}}\cdot \mathbf{u}]]=0,$$

$$\frac{\text{d}\mathbf{x}}{\text{d}t}=(\hat{\mathbf{n}}\cdot \mathbf{u})\hat{\mathbf{n}},$$

$$\text{mbox}Stable:M=\frac{{\mu }^{-}}{{\mu }^{+}}>1.$$

$$\omega =\frac{k}{L}\left(U\frac{{\mu }^{+}-{\mu }^{-}}{{\mu }^{+}+{\mu }^{-}}-\frac{\gamma b^2{k}^2}{12L^2\left({\mu }^{+}+{\mu }^{-}\right)}\right),$$

$${\mathbf{F}}^{\pm }=-{\mathbb{K}}^{\pm }\nabla {\mathbf{\Phi }}^{\pm },\nabla \cdot {\mathbf{F}}^{\pm }=0,$$

$$\mathbb{K}=\left(\begin{array}[]{cc}K_{h}&K_{eo}\\ K_{eo}&K_{e}\end{array}\right).$$

$$-\nabla \mathbf{\Phi }\cdot \mathbf{F}=\nabla {\mathbf{\Phi }}^{\mathsf{T}}\mathbb{K}\nabla \mathbf{\Phi }>0.$$

$$[[\mathbf{\Phi }]]=(\gamma \kappa ,0{)}^{\mathsf{T}},[[\hat{\mathbf{n}}\cdot \mathbf{F}]]=0.$$

$$\omega =\frac{k}{L}\left(F-\gamma G\frac{k^2}{L^2}\right),$$

$$F=U\frac{[[K_{eo}]]\{K_{eo}\}-[[K_h]]\{K_e\}}{\det \{\mathbb{K}\}}+2I\frac{K_h^{+}{K}_{eo}^{-}-K_h^{-}{K}_{eo}^{+}}{\det \{\mathbb{K}\}},$$

$$G=\frac{K_h^{+}\det {\mathbb{K}}^{-}+K_h^{-}\det {\mathbb{K}}^{+}}{\det \{\mathbb{K}\}},$$

$${\lambda }_{cr}=\pi b\sqrt{\frac{\gamma }{3U[[\mu ]]+6I[[\epsilon \zeta ]]/\{\sigma \}}}.$$

$$\text{mbox}Stable:\frac{M-1}{M+1}>A\frac{SZ-1}{SZ+1},$$

$$S=\frac{{\epsilon }^{-}}{{\epsilon }^{+}},Z=\frac{{\zeta }^{-}}{{\zeta }^{+}},M=\frac{{\mu }^{-}}{{\mu }^{+}},A=\frac{I(-\overline{\epsilon \zeta })}{U\bar{\mu }\bar{\sigma }},$$

$$I_{cr}\approx \frac{U}{2}\frac{{\mu }^o}{{\mu }^w}\frac{{\sigma }^w}{K_{eo}^w}.$$

$$E_{cr}\approx \frac{U}{4}\frac{{\mu }^o}{{\mu }^w}\frac{{\sigma }^w}{{\sigma }^o}\frac{1}{K_{eo}^w}.$$

$$\omega =-\frac{U}{r}+\frac{k}{r}\left(F-\gamma G\frac{(k^2-1)}{r^2}\right),$$

Full text

Interfacial instability is the precursor to pattern formation in a variety of physical and chemical processes Kessler et al. (1988); Pelcé and Libchaber (2012). This fascinating topic covers a broad range of phenomena such as dendritic growth due to the Mullins-Sekerka instability in solidification Mullins and Sekerka (1963); *Mullins1964; Langer (1980), fractal growth due to diffusion-limited aggregation Witten and Sander (1983) or metal electrodeposition Brady and Ball (1984) in fluid flows Bazant et al. (2003), crease formation and wrinkling of combustion fronts due to the Darrieus-Landau instability Darrieus ; *landau1944theory; Sivashinsky (1983); *matalon2007intrinsic, and viscous fingering in Hele-Shaw cells Bensimon et al. (1986) and porous media Homsy (1987) due to the Saffman-Taylor instability Saffman and Taylor (1958); Chuoke et al. (1959).

Interfacial instabilities are usually undesirable, but difficult to control. In secondary oil recovery, viscous fingering of injected liquids leads to nonuniform displacement and residual trapping of oil Homsy (1987); Gorell and Homsy (1983), and dendritic growth is a major safety concern for metal anodes in rechargeable batteries Xu et al. (2014). There are signs, however, that instability may be avoided if the interface is driven by multiple opposing forces. For instance, it was recently observed that dendritic growth can be suppressed in charged porous media Han et al. (2016) if preceded by deionization shock wave Mani and Bazant (2011), whose stable propagation in cross flow also enables water purification by shock electrodialysis Deng et al. (2013); Schlumpberger et al. (2015).

Here, we consider the interfacial stability of two immiscible electrolytes in a Hele-Shaw cell where the interface is set into motion by both the pressure-driven and electro-osmotic flows. Remarkably, we find that electrokinetic coupling influence interfacial stability and, under certain conditions, can eliminate viscous fingering. This phenomenon illustrates the rich physics of “Vector Laplacian Growth” (VLG), a general mathematical model of interfacial dynamics driven by the gradient of a vector-valued potential function through a generalized mobility tensor. The “one-sided” VLG model (with field gradients only on one side of the interface) is known to be unstable, leading to fractal patterns, during growth Bazant et al. (2003) and stable, resulting in smooth collapse, during retreat Bazant (2006); Rycroft and Bazant (2016). Our theory shows that stable growth is also possible, if field gradients exist on both sides of the interface.

In the classical viscous fingering problem, the fluid flow in a Hele-Shaw cell can be approximated as quasi two-dimensional if the cell gap, $b$, is much smaller than the lateral dimension, $L$ (see fig 1). In this case, the gap-averaged velocity of each fluid is given by:

[TABLE]

where ‘$-$’ and ‘$+$’ superscripts denote invading and receding fluids, $\mu$ and $p$ are viscosity and pressure of each fluid, and $\nabla$ is the in-plane gradient operator. At the interface, the pressure jump is given by the Young-Laplace equation, while the normal velocity is continuous:

[TABLE]

where $\llbracket a\rrbracket\equiv a^{+}-a^{-}$ denotes the jump of variable ‘$a$’ across the interface, $\gamma$ is the surface tension, and $\kappa$ is the in-plane curvature. More generally, these conditions must be modified to take the finite lubrication film thickness into account if the receding fluid is perfectly wetting Park and Homsy (1984); Homsy (1987). The interface moves with the local fluid velocity,

[TABLE]

and the far-field flow is uniform.

Linear stability analysis of Eqs. (1)–(3), initiated by Chuoke et al. Chuoke et al. (1959) and Saffman and Taylor Saffman and Taylor (1958), reveals that stable displacement is only possible if the advancing fluid is more viscous:

[TABLE]

In the opposite case, $M<1$, the interface is unstable to perturbations of sufficiently long wavelength, and the less viscous fluid forms “fingers” of lower resistance through the more viscous fluid. Specifically, the growth rate, $\omega$, of a normal mode $\delta y\sim\exp(ikx/L+\omega t)$ satisfies the dispersion relation Homsy (1987); Bensimon et al. (1986):

[TABLE]

where $k$ is the wavenumber. Perturbation wavelengths longer than $\lambda_{cr}=\pi b\sqrt{\gamma/3U(\mu^{+}-\mu^{-})}$ are unstable, and the maximum growth rate arises for $\lambda_{m}=\sqrt{3}\lambda_{cr}$.

Most materials naturally acquire charge in aqueous solutions from the dissociation of surface groups, such as silanol Kirby and Hasselbrink (2004); Hunter (1988), for glass in Hele-Shaw cells or silicate minerals in underground reservoirs. The screening of surface charge by mobile ions leads to the formation of electric double layers (EDL) and associated electrokinetic phenomena Russel et al. (1989); *lyklema1995fundamentals. An electric field parallel to the charged surface acts on EDL charge to drive “electro-osmotic” flow $\mathbf{u}_{eo}$, while pressure-driven flow drives “streaming current” $\mathbf{i}_{sc}$ due to the advection of EDL charge (Fig. 1). For typical situations of fixed surface charge, the electrokinetic response is linear in the driving forces, i.e. $\mathbf{u}_{eo}=-K_{eo}\nabla\phi$ and $\mathbf{i}_{sc}=-K_{sc}\nabla p$, where $\phi$ is the electrostatic potential. The electro-osmotic mobility, $K_{eo}$, and the streaming conductance, $K_{sc}$, satisfy Onsager’s reciprocal relation Onsager (1931); Bazant (2016), $K_{eo}=K_{sc}$, and for thin EDL (gaps, $b\sim 0.1-1\>\text{mm}$, much larger than the EDL thickness, $\lambda_{D}\sim 1-10\>\text{nm}$), are given by the Helmholtz-Smoluchowski relation Lyklema (1995), $K_{eo}=-\varepsilon\zeta/\mu$, where $\varepsilon$ is the electrolyte permittivity and $\zeta$ is the potential difference across the EDL Lyklema (1995); Russel et al. (1989); Hunter (1988).

When linear electrokinetic phenomena are considered, a VLG model can thus be written in terms of a tensorial flux, $\mathbf{F}=(\mathbf{u},\mathbf{i})^{\mathsf{T}}$, proportional to the gradient of a vector-valued potential, $\mathbf{\Phi}=(p,\phi)^{\mathsf{T}}$:

[TABLE]

where $\mathbb{K}$ is the electrokinetic mobility tensor:

[TABLE]

$K_{h}=b^{2}/12\mu$ is the hydraulic Darcy conductivity, and $K_{e}=\sigma$ is the electrical Ohmic conductivity of the cell. The Second Law of Thermodynamics requires positive definite $\mathbb{K}$ to ensure positive dissipation rate De Groot and Mazur (2013); Peters et al. (2016), i.e.:

[TABLE]

At the interface, the pressure and total velocity satisfy the jump conditions given by Eq. (2), while the potential and normal component of the total current are continuous, which can be compactly expressed as:

[TABLE]

Far from the interface in a planar geometry, the fluxes are assumed to be uniform, $\lim_{y\rightarrow\pm\infty}\mathbf{F}_{y}=\mathbf{F}_{\infty}=(U,I)^{\mathsf{T}}$. Equations (6) and (9), along with the kinematic condition (3), determine the interface motion.

As for classical problem, we consider the linear stability of a planar interface subjected to a sinusoidal perturbation, $\delta y\sim\exp(ikx/L+\omega t)$, and seek solutions of the form $\mathbf{\Phi}^{\pm}=\mathbf{\Phi}^{\pm}_{0}+\epsilon\mathbf{\Phi}^{\pm}_{1}$ in the limit of $\epsilon\ll 1$. From Eq. (6), the base state is linear, i.e. $\mathbf{\Phi}^{\pm}_{0}=-\mathbb{K}^{\pm^{-1}}\mathbf{F}_{\infty}(y-Ut)$, while $\mathbf{\Phi}^{\pm}_{1}=\mathbf{A}^{\pm}_{1}\exp(ikx/L+\omega t)\exp(\mp k(y-Ut)/L)$, where $\mathbf{A}^{\pm}_{1}$ are evaluated using the jump condition (9). Applying the kinematic condition (3) then yields the growth rate:

[TABLE]

where $F$ and $G$ are given by

[TABLE]

and $\{a\}\equiv a^{+}+a^{-}$. Note that the classical dispersion relation (5) is recovered in the absence of electrokinetic phenomena, $K_{eo}^{\pm}=0$. From the Second Law (8), it follows that $G>0$, ensuring that surface tension effects are stabilizing. Therefore, $F<0$ is a sufficient condition for stability. For $F>0$, a perturbation of wavelength longer than $\lambda_{cr}=2\pi\sqrt{\gamma G/F}$ is unstable, and $\lambda_{m}=\sqrt{3}\lambda_{cr}$ is the most unstable wavelength.

To simplify equation (11), we note that the electrokinetic coupling coefficient Van Der Heyden et al. (2006); *VanDerHeyden2007, $\alpha=K_{eo}^{2}/K_{h}K_{e}$, is typically small, while $0\leq\alpha<1$ from Eq. (8). For $\alpha\ll 1$, the critical wavelength may be approximated as:

[TABLE]

While the classical instability is controlled by the viscosity ratio (Eq. (4)), our theory predicts that the injection ratio, $I/U$, can be tuned independently to control interfacial stability (Fig. 2a):

[TABLE]

in terms of the following dimensionless ratios:

[TABLE]

where the over-bar indicates average values, e.g. $\overline{\varepsilon\zeta}=\left(\varepsilon^{+}\zeta^{+}+\varepsilon^{-}\zeta^{-}\right)/2$. Electrokinetic effects require $SZ\neq 1$, and stability is possible if the injection ratio is larger than $\left|(M-1)(SZ+1)/(SZ-1)(M+1)\right|$, and has the “correct” sign, depending on the magnitude of $SZ$ (see Fig. 2). Above a critical injection ratio, $A_{cr}=\left|(SZ+1)/(SZ-1)\right|$, stability is entirely determined by the sign of injection ratio and is, remarkably, independent of the viscosity ratio. Physically, negative injection ratios denote opposite direction of current and flow. These observations are illustrated in figure 2b–d.

Motivated by secondary oil recovery, it is interesting to consider the limit when $M\ll 1$, $SZ\gg 1$, and $R=\sigma^{-}/\sigma^{+}\gg 1$, e.g. when water is pushing oil toward extraction wells. In this case, negative current injection shifts the critical wavelength to longer values and reduces viscous fingering. Stable displacement is possible if $I>I_{cr}$, where:

[TABLE]

For $U=1\>\text{mm}/\text{min}$ in a $1\>\text{mM}$ KCl solution with $\zeta=-50\>\text{mV}$ and $\mu^{w}/\mu^{o}=0.1$, the critical current is fairly small, $I_{cr}\sim 4\>\text{mA}/\text{cm}^{2}$, but a large critical electric field is required to drive this current across the poorly conducting oil region:

[TABLE]

Even for a modest value of $\sigma^{w}/\sigma^{o}=10$, complete suppression of viscous fingering requires $E_{cr}\sim 150\>\text{V}/\text{cm}$. The required voltage could be lowered by reducing the conductivity ratio or electrode separations. Nonetheless, partial stabilization (with enhanced oil recovery) is still viable with electric fields below the critical value.

To further support our theory, we numerically solve the VLG model using the Voronoi Interface Method Guittet et al. (2015) to discretize the conservation equations (6) subjected to the interface jump conditions (9) while utilizing the level-set framework Osher and Sethian (1988); *osher2006level; *sethian1999level to represent the moving interface. Furthermore, we use dynamically adaptive quadtree grids Min and Gibou (2007) as well as parallel algorithms Mirzadeh et al. (2016) for fast and high-fidelity simulations. Figure 3 illustrates the interfacial dynamics of an initial perturbation for unfavorable, 3-(a,b), and favorable, 3-(c,d), viscosity ratios. As predicted, interfacial stability can be manipulated by adjusting the injection ratio.

It is straight-forward to extend the analysis for the radial Hele-Shaw cell geometry Paterson (1981); Homsy (1987) where an invading fluid is injected at a point to push the second fluid outward. If the interface is initially assumed to be circular, the growth rate of an azimuthal perturbation of the form $\delta r\sim\exp\left(ik\theta+\omega t\right)$, is given via:

[TABLE]

where $r$ is the initial radius and $F$ and $G$ are still given by equation (11). Once again, $F<0$ is a sufficient condition for stability. Therefore, the stability estimate in equation (13) could be used in radial geometry if velocity, $U$, and current density, $I$, are replaced by the total flow rate, $Q_{0}$, and total current, $I_{0}$, respectively. Figure 4 illustrates numerical simulation of interface evolution in a radial Hele-Shaw cell geometry for an unstable viscosity ratio of $M=0.01$. The instability is entirely suppressed when current is injected in the opposite direction.

The possibility of manipulating interfacial instabilities is quite exciting. The idea of controlling viscous fingering using cell geometry has been recently discussed Al-Housseiny et al. (2012); *Al-Housseiny2013a, and our framework introduces many other degrees of freedom, such as the placement of electrodes, dielectric or conducting boundaries, and surface coatings or gate voltages to modify local zeta potentials. For a given geometry, dynamical control of fingering instabilities may also be possible, by adjusting potentials and pressures with real-time feedback from currents and flow rates.

Hele-Shaw cell experiments could be used to check these predictions and test the validity of our assumptions. Since electrokinetic phenomena depend on which liquid is in contact with the surface, it may be necessary to extend the model for lubrication films and gravity currents Eames et al. (2005), which would require more complicated depth-averaging and electrokinetics at the liquid-liquid interface Pascall and Squires (2011), perhaps amenable to conformal-map dynamics Bazant et al. (2003); Robinson and Eames (2017). We have also neglected non-linear electro-hydrodynamic effects Melcher and Taylor (1969); *Melcher1981; *melcher1963field which might cause interfacial instabilities at higher electric fields in large channels Taylor and McEwan (1965); Lin et al. (2004); *Chen2005; Druzgalski et al. (2013). We also assume finite electrical resistance in each phase, which could exclude traditional liquid pairs, such as water/silicon oil and air/glycerol, although some poorly conducting regions may have sufficient ionic or electronic conductivity to pass at least a transient current consistent with the model. The model could be extended to include interfacial capacitance, and the resulting “RC time” for charge accumulation might be longer than the instability growth time, especially for large resistive domains.

Our model is directly applicable to interfaces between two immiscible electrolyte solutions (ITIES) that support charge-transfer reactions Girault (1987). Examples include aqueous electrolytes, e.g. LiCl, in contact with solutions of lipophilic salts in organic solvents, e.g. TBATPB in nitrobenzene Senda et al. (1991); *girault1993charge; *Samec2004; *girault2010electrochemistry. Recent interest in ITIES was spurred by electro-wetting Monroe et al. (2006) for electro-variable optics Edel et al. (2016), but tunable fingering under confinement could lead to different applications.

Although our theory is for immiscible electrolytes, it may also describe diffuse interfaces involving strong ion concentration gradients, e.g. deionization shocks in charged porous media Mani and Bazant (2011)or pH fronts in electrokinetic remediation of contaminated soil Probstein and Hicks (1993); *Shapiro1993. Since the $\zeta$-potential is a function of pH and salt concentration Kirby and Hasselbrink (2004); Hunter (1988), it may also be possible to observe some of the stabilizing effects with miscible solutions, perhaps in charged porous media such as glass frit or Hele-Shaw cells packed with silica beads. Finally, we caution that viscous fingering is more complicated in porous media than in Hele-Shaw cells, due to permeability variations, capillary effects, and surface wettability Lenormand et al. (1988); Zhao et al. (2016). Since electrokinetic couplings derive from surfaces, we expect strong dependence on surface wettability whereas permeability variations might have limited impact due to disproportionate scaling of hydraulic and electro-osmotic mobilities with the pore size, possibly resulting in a more uniform displacement. Nonetheless, further investigation is required to quantify the degree to which electrokinetic phenomena can control interfacial stability in porous media.

This work was supported by a seed grant from the MIT Energy Initiative with computational resources from the Texas Advanced Computing Center (TACC) at The University of Texas at Austin and the Extreme Science and Engineering Discovery Environment (XSEDE), supported by grant ACI-1548562 from the National Science Foundation. The authors thank Amir A. Pahlavan and Charles W. Monroe for useful discussions.