2.2.1 Capturing salinity variations

In order to study the effect of varying the concentrations of monovalent and divalent ions simultaneously, since this is typical for water-flood experiments, we introduce the parameter $\unicode[STIX]{x1D70F}$ such that

(2.30) $$\begin{eqnarray}\displaystyle c_{\infty }^{1,2}=c_{low}^{1,2}+10^{\unicode[STIX]{x1D70F}}(c_{high}^{1,2}-c_{low}^{1,2}), & & \displaystyle\end{eqnarray}$$

where $c_{low}^{1,2}$ are the concentrations of monovalent and divalent ions, respectively, in the low-salinity regime, and $c_{high}^{1,2}$ are the concentrations of monovalent and divalent ions, respectively, in the high-salinity regime. By varying $\unicode[STIX]{x1D70F}$ from $-\infty$ to 0, we can study the model continuously between the low-salinity and high-salinity regimes using a single parameter. We choose a power law for $\unicode[STIX]{x1D70F}$ , as typically the salinities span multiple orders of magnitude.

3.3.1 Regions I and III

In regions I and III, we expand $\unicode[STIX]{x1D719}$ and $d^{\pm ,2\pm }$ in powers of $\unicode[STIX]{x1D6FF}$ to find that the leading-order contributions, $\unicode[STIX]{x1D719}_{0}$ and $d_{0}^{\pm ,2\pm }$ , satisfy

(3.9) $$\begin{eqnarray}\displaystyle d_{0}^{-}\text{e}^{\unicode[STIX]{x1D719}_{0}}-d_{0}^{+}\text{e}^{-\unicode[STIX]{x1D719}_{0}}+2c^{\ast }(d_{0}^{2-}\text{e}^{2\unicode[STIX]{x1D719}_{0}}-d_{0}^{2+}\text{e}^{-2\unicode[STIX]{x1D719}_{0}})=0. & & \displaystyle\end{eqnarray}$$

Since (3.9) is an algebraic expression relating $\unicode[STIX]{x1D719}_{0}$ and $d_{0}^{\pm ,2\pm }(x,t)$ , it follows that $\unicode[STIX]{x1D719}_{0}=\unicode[STIX]{x1D719}_{0}(x,t)$ , and therefore at leading order (3.3)–(3.5) become

(3.10) $$\begin{eqnarray}\displaystyle & \displaystyle Pe^{\pm }\left[(d_{0}^{\pm }\text{e}^{\mp \unicode[STIX]{x1D719}_{0}})_{t}+(d_{0}^{\pm }\text{e}^{\mp \unicode[STIX]{x1D719}_{0}})_{x}\int _{0}^{1}u_{w}\,\text{d}z\right]=(d_{0_{x}}^{\pm }\text{e}^{\mp \unicode[STIX]{x1D719}_{0}})_{x}, & \displaystyle\end{eqnarray}$$

(3.11) $$\begin{eqnarray}\displaystyle & \displaystyle Pe^{2\pm }\left[(d_{0}^{2\pm }\text{e}^{\mp 2\unicode[STIX]{x1D719}_{0}})_{t}+(d_{0}^{2\pm }\text{e}^{\mp 2\unicode[STIX]{x1D719}_{0}})_{x}\int _{0}^{1}u_{w}\,\text{d}z\right]=(d_{0_{x}}^{2\pm }\text{e}^{\mp 2\unicode[STIX]{x1D719}_{0}})_{x}, & \displaystyle\end{eqnarray}$$

for $0<x<x_{r}^{\ast }$ and $x_{f}^{\ast }<x<1$ , where $x_{r}^{\ast }=x_{r}/L$ and $x_{f}^{\ast }=x_{f}/L$ . Equations (3.10) and (3.11) are subject to $d_{0}^{\pm ,2\pm }(1,t)=d_{0}^{\pm ,2\pm }(0,t)=1+10^{\unicode[STIX]{x1D70F}}(\bar{c}^{1,2}-1)$ and

(3.12) $$\begin{eqnarray}\displaystyle \frac{(d_{0}^{+})_{x}\,\text{e}^{-\unicode[STIX]{x1D719}_{0}}}{Pe^{+}}-\frac{(d_{0}^{-})_{x}\,\text{e}^{\unicode[STIX]{x1D719}_{0}}}{Pe^{-}}+2c^{\ast }\left[\frac{(d_{0}^{2+})_{x}\,\text{e}^{-2\unicode[STIX]{x1D719}_{0}}}{Pe^{2+}}-\frac{(d_{0}^{2-})_{x}\,\text{e}^{2\unicode[STIX]{x1D719}_{0}}}{Pe^{2-}}\right]=0, & & \displaystyle\end{eqnarray}$$

at $x=0$ , from (2.57)–(2.59). We will find that the relevant solutions are $\unicode[STIX]{x1D719}_{0}$ constant, which we take to be zero without loss of generality, and $d_{0}^{\pm ,2\pm }\equiv 1+10^{\unicode[STIX]{x1D70F}}(\bar{c}^{1,2}-1)$ by matching with the solutions in regions IV and V.

However, after imposing (2.60) and solving $\unicode[STIX]{x1D719}_{zz}=0$ in the clay phase, we find that $\unicode[STIX]{x1D719}_{z}$ becomes unbounded as $z\rightarrow -\infty$ , which contradicts (2.63). Hence, there is a boundary layer at the capillary wall of size $\unicode[STIX]{x1D716}_{2}$ (to be determined). We rescale $z=\unicode[STIX]{x1D716}_{2}Z$ . From (2.60) and (2.61), we observe that the rescaling $\unicode[STIX]{x1D719}=\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1}\unicode[STIX]{x1D6F7}$ is required to balance the terms in the boundary conditions.

We expand $d^{\pm ,2\pm }\sim \tilde{d_{0}}^{\pm ,2\pm }+\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1}\tilde{d_{1}}^{\pm ,2\pm }+O(\unicode[STIX]{x1D716}_{2}^{2}\unicode[STIX]{x1D709}^{-2})$ (assuming that $\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1}\ll 1$ ) and substitute into (3.2) to get

(3.13) $$\begin{eqnarray}\displaystyle \frac{\unicode[STIX]{x1D6FF}\unicode[STIX]{x1D709}}{\unicode[STIX]{x1D716}_{2}^{2}}\unicode[STIX]{x1D6F7}_{ZZ}=[d_{0}^{+}+d_{0}^{-}+4c^{\ast }(d_{0}^{2+}+d_{0}^{2-})]\unicode[STIX]{x1D6F7}+\tilde{d_{1}}^{-}-\tilde{d_{1}}^{+}+2c^{\ast }(\tilde{d_{1}}^{2-}-\tilde{d_{1}}^{2+})+O\left(\frac{\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}\right), & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

where we have used the fact that $\tilde{d_{0}}^{\pm ,2\pm }=d_{0}^{\pm ,2\pm }$ by matching with the solution for $d^{\pm ,2\pm }$ in the outer region. To balance the leading-order terms in (3.13), we substitute $d_{0}^{-,2-}=d_{0}^{+,2+}$ and set

(3.14) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D716}_{2}=\left(\frac{\unicode[STIX]{x1D6FF}\unicode[STIX]{x1D709}}{2(d_{0}^{+}+4c^{\ast }d_{0}^{2+})}\right)^{1/2}=\frac{\unicode[STIX]{x1D706}}{R}, & & \displaystyle\end{eqnarray}$$

where

(3.15) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D706}=\left(\frac{\unicode[STIX]{x1D700}_{w}k_{B}T}{2(c_{\infty }^{1}+4c_{\infty }^{2})q^{2}}\right)^{1/2} & & \displaystyle\end{eqnarray}$$

is the Debye length. From (3.14) we see that $\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1}=O(\unicode[STIX]{x1D6FF}^{1/2}\unicode[STIX]{x1D709}^{-1/2})$ . Since $\unicode[STIX]{x1D6FF}\ll \unicode[STIX]{x1D709}$ , the assumption we made that $\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1}\ll 1$ is valid, so our asymptotic analysis is consistent. The dimensionless parameter $\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1}$ is the ratio of the adsorption depth to the Debye length.

Since $\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1}\gg \unicode[STIX]{x1D6FF}$ , and the next non-zero terms in the expansions for $d^{\pm ,2\pm }$ in the outer region away from the capillary wall are $O(\unicode[STIX]{x1D6FF})$ , it follows that the $O(\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1})$ terms, $\tilde{d_{1}}^{\pm ,2\pm }$ , in the expansion for $d^{\pm ,2\pm }$ near the capillary wall are zero. Hence, with this definition of $\unicode[STIX]{x1D716}_{2}$ , equation (3.13) reads

(3.16) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6F7}_{ZZ}=\unicode[STIX]{x1D6F7}, & & \displaystyle\end{eqnarray}$$

for $0<x<x_{r}^{\ast }$ , $0<z<\infty$ , and $x_{f}^{\ast }<x<1$ , $0<z<\infty$ at leading order in $\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1}$ , subject to

(3.17) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6F7}_{Z}\rightarrow 0,\quad \text{as }Z\rightarrow \infty , & & \displaystyle\end{eqnarray}$$

(3.18a,b ) $$\begin{eqnarray}\displaystyle [\unicode[STIX]{x1D6F7}]_{Z=0^{-}}^{Z=0^{+}}=0,\quad \unicode[STIX]{x1D6F7}_{Z}=2{s_{c}}_{0}^{+}+{s_{c}}_{0}-1,\quad \text{at }Z=0. & & \displaystyle\end{eqnarray}$$

Equation (3.17) is the matching condition required to ensure that the potential in the boundary layer, $\unicode[STIX]{x1D6F7}$ , matches with the outer solution, $\unicode[STIX]{x1D719}_{0}=0$ . Note that we have used the fact that $\unicode[STIX]{x1D719}_{z}=0$ in the clay phase to derive (3.18). We solve (3.16) subject to (3.17) and (3.18) to get

(3.19) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6F7}=-(2{s_{c}}_{0}^{+}+{s_{c}}_{0}-1)\text{e}^{-Z}. & & \displaystyle\end{eqnarray}$$

From (3.19), we see that there is a boundary layer of the order of the Debye length across which the electric potential decays exponentially due to screening from counter-ions in the water phase.

3.3.2 Region II

We now turn our attention to finding the potential in region II. Since we expect the thickness of the wetting film to be partly determined by the electrostatic forces occurring between the two charged surfaces, we suppose that $h^{\ast }=O(\unicode[STIX]{x1D716}_{2})$ . (This is the canonical scaling. If $h^{\ast }$ is larger than $O(\unicode[STIX]{x1D716}_{2})$ , there is a bulk region between the surfaces in which $\unicode[STIX]{x1D719}=0$ , and therefore there is no electrical interaction between them.) Therefore, we rescale $h=\unicode[STIX]{x1D716}_{2}H$ . To balance the terms in (2.52) and (2.61), we rescale $u_{w}=\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}U_{w}$ , $w_{w}=\unicode[STIX]{x1D716}_{2}^{2}\unicode[STIX]{x1D709}W_{w}$ and $t=\unicode[STIX]{x1D716}_{2}^{-1}\unicode[STIX]{x1D709}^{-1}T$ . In terms of $Z$ , $H$ , $\unicode[STIX]{x1D6F7}$ , $U_{w}$ , $W_{w}$ and $T$ , equations (3.2) and (3.6)–(3.8) read

(3.20) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6F7}_{ZZ}=\frac{\unicode[STIX]{x1D709}}{\unicode[STIX]{x1D716}_{2}}\frac{d^{-}-d^{+}+2c^{\ast }(d^{2-}-d^{2+})}{2(d_{0}^{+}+4c^{\ast }d_{0}^{2+})}+\frac{d^{+}+d^{-}+4c^{\ast }(d^{2+}+d^{2-})}{2(d_{0}^{+}+4c^{\ast }d_{0}^{2+})}\unicode[STIX]{x1D6F7}+O\left(\frac{\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}\right), & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

(3.21) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D709}^{2}\frac{\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}Pe^{+}\left[Hd_{T}^{+}+d_{x}^{+}\int _{0}^{H}U_{w}\,\text{d}Z\right] & = & \displaystyle \left[d_{x}^{+}\int _{0}^{H}\left(1-\frac{\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}\unicode[STIX]{x1D6F7}\right)\text{d}Z\right]_{x}\nonumber\\ \displaystyle & & \displaystyle -\,\frac{\unicode[STIX]{x1D6FF}\unicode[STIX]{x1D709}\,Pe^{2+}}{c^{\ast }}\left(s^{\ast }\frac{\text{D}s_{o}^{+}}{\text{D}T}+\frac{\unicode[STIX]{x2202}s_{c}^{+}}{\unicode[STIX]{x2202}T}\right)+O\left(\left(\frac{\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}\right)^{2}\right),\nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

(3.22) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D709}^{2}\frac{\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}Pe^{2+}\left[Hd_{T}^{2+}+d_{x}^{2+}\int _{0}^{H}U_{w}\,\text{d}Z\right] & = & \displaystyle \left[d_{x}^{2+}\int _{0}^{H}\left(1-\frac{2\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}\unicode[STIX]{x1D6F7}\right)\text{d}Z\right]_{x}\nonumber\\ \displaystyle & & \displaystyle -\,\frac{\unicode[STIX]{x1D6FF}\unicode[STIX]{x1D709}\,Pe^{2+}}{c^{\ast }}\left(s^{\ast }\frac{\text{D}s_{o}^{+}}{\text{D}T}+\frac{\unicode[STIX]{x2202}s_{c}^{+}}{\unicode[STIX]{x2202}T}\right)+O\left(\left(\frac{\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}\right)^{2}\right),\nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

(3.23) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D709}^{2}\frac{\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}Pe^{k-}\left[Hd_{T}^{k-}+d_{x}^{k-}\int _{0}^{H}U_{w}\,\text{d}Z\right]=\left[d_{x}^{k-}\int _{0}^{H}\left(1+\frac{k\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}\unicode[STIX]{x1D6F7}\right)\text{d}Z\right]_{x}+O\left(\left(\frac{\unicode[STIX]{x1D716}_{2}}{\unicode[STIX]{x1D709}}\right)^{2}\right), & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

for $k=1,2$ . We expand $\unicode[STIX]{x1D6F7}$ and $d^{\pm ,2\pm }$ in powers of $\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1}$ to find that, at leading order, equations (3.20)–(3.23) read

(3.24) $$\begin{eqnarray}\displaystyle d_{0}^{-}-d_{0}^{+}+2c^{\ast }(d_{0}^{2-}-d_{0}^{2+})=0, & & \displaystyle\end{eqnarray}$$

(3.25a,b ) $$\begin{eqnarray}\displaystyle ({d_{0}^{\pm }}_{x}H)_{x}=0,\quad ({d_{0}^{2\pm }}_{x}H)_{x}=0. & & \displaystyle\end{eqnarray}$$

We have used the fact that $\unicode[STIX]{x1D6FF}\unicode[STIX]{x1D709}=O((\unicode[STIX]{x1D716}_{2}/\unicode[STIX]{x1D709})^{2}\unicode[STIX]{x1D709}^{2})\ll (\unicode[STIX]{x1D716}_{2}/\unicode[STIX]{x1D709})^{2}$ , so that the terms involving the surface concentrations in (3.21) and (3.22) are neglected at leading order. We will find that the relevant solutions to (3.24) and (3.25) are the constant solutions $d_{0}^{\pm ,2\pm }\equiv 1+10^{\unicode[STIX]{x1D70F}}(\bar{c}^{1,2}-1)$ , by matching with the solutions in regions IV and V.

At the next order in $\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1}$ , equation (3.20) reads

(3.26) $$\begin{eqnarray}\displaystyle {\unicode[STIX]{x1D6F7}_{0}}_{ZZ}=\unicode[STIX]{x1D6F7}_{0}+\frac{d_{1}^{-}-d_{1}^{+}+2c^{\ast }(d_{1}^{2-}-d_{1}^{2+})}{2(d_{0}^{+}+4c^{\ast }d_{0}^{2+})}. & & \displaystyle\end{eqnarray}$$

We solve this subject to the boundary conditions

(3.27a,b ) $$\begin{eqnarray}\displaystyle [\unicode[STIX]{x1D6F7}_{0}]_{Z=0^{-}}^{Z=0^{+}}=0,\quad {\unicode[STIX]{x1D6F7}_{0}}_{Z}=2{s_{c}}_{0}^{+}+{s_{c}}_{0}-1,\quad \text{at }Z=0, & & \displaystyle\end{eqnarray}$$

(3.28a,b ) $$\begin{eqnarray}\displaystyle [\unicode[STIX]{x1D6F7}_{0}]_{Z=H^{-}}^{Z=H^{+}}=0,\quad {\unicode[STIX]{x1D6F7}_{0}}_{Z}|^{Z=H^{-}}=-s^{\ast }(2{s_{o}}_{0}^{+}+{s_{o}}_{0}-1),\quad \text{at }Z=H, & & \displaystyle\end{eqnarray}$$

to get

(3.29) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6F7}_{0} & = & \displaystyle -\frac{d_{1}^{-}-d_{1}^{+}+2c^{\ast }(d_{1}^{2-}-d_{1}^{2+})}{2(d_{0}^{+}+4c^{\ast }d_{0}^{2+})}+(2{s_{c}}_{0}^{+}+{s_{c}}_{0}-1)\sinh Z\nonumber\\ \displaystyle & & \displaystyle -\,\frac{(2{s_{c}}_{0}^{+}+{s_{c}}_{0}-1)\cosh H+s^{\ast }(2{s_{o}}_{0}^{+}+{s_{o}}_{0}-1)}{\sinh H}\cosh Z.\end{eqnarray}$$

We have used the fact that $\unicode[STIX]{x1D719}_{z}=0$ in the oil and clay phases to derive the boundary conditions (3.27) and (3.28). We have also used the fact that $H_{x}$ is small in order to make the approximation $\boldsymbol{n}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D6F7}\sim -\unicode[STIX]{x1D6F7}_{Z}$ on the oil surface, which is required to derive (3.28).

The first-order terms in (3.21)–(3.23) read

(3.30a,b ) $$\begin{eqnarray}\displaystyle [(d_{1}^{\pm })_{x}H]_{x}=0,\quad [(d_{1}^{2\pm })_{x}H]_{x}=0. & & \displaystyle\end{eqnarray}$$

In order to verify the solution $d_{0}^{\pm ,2\pm }\equiv 1+10^{\unicode[STIX]{x1D70F}}(\bar{c}^{1,2}-1)$ throughout the water phase, and $d_{1}^{\pm ,2\pm }\equiv 0$ in the thin film, we must match together the solutions in regions I, II and III via regions IV and V. A schematic diagram of the matching process is shown in figure 4. Towards the left of region IV, the oil and clay surfaces are far enough apart for the bulk water solution, $\unicode[STIX]{x1D719}_{0}\equiv 0$ , $d_{0}^{\pm ,2\pm }\equiv 1+10^{\unicode[STIX]{x1D70F}}(\bar{c}^{1,2}-1)$ , to be valid. We assume that, in region IV, the meniscus is approximately circular, which, dimensionally, reads $h\sim (x-x_{c})^{2}/2R$ . In terms of the dimensionless variables ( $x=L\hat{x}$ , $h=\unicode[STIX]{x1D716}_{2}R{\hat{H}}$ , and subsequently dropping the hat notation for convenience) this reads $H\sim (1/2)\unicode[STIX]{x1D716}_{1}^{-2}\unicode[STIX]{x1D716}_{2}^{-1}(x-x_{c}^{\ast })^{2}$ for $|x-x_{c}^{\ast }|\ll 1$ , where $x_{c}^{\ast }=x_{c}/L$ is the dimensionless apparent contact point. Thus, for $|x-x_{c}^{\ast }|\gg \unicode[STIX]{x1D716}_{1}\unicode[STIX]{x1D716}_{2}^{1/2}$ , we find that $H\gg 1$ and therefore the relevant solution is $\unicode[STIX]{x1D719}_{0}\equiv 0$ and $d_{0}^{\pm ,2\pm }\equiv 1+10^{\unicode[STIX]{x1D70F}}(\bar{c}^{1,2}-1)$ , with a boundary layer on the clay surface in which $\unicode[STIX]{x1D6F7}$ is given by (3.19).

Towards the right of region IV, the gradient of the oil–water interface is small so (3.29) applies. This approximation is valid provided $H_{x}\ll \unicode[STIX]{x1D716}_{1}^{-1}\unicode[STIX]{x1D716}_{2}^{-1}$ , which holds for $|x-x_{c}^{\ast }|\ll \unicode[STIX]{x1D716}_{1}$ , assuming the circular shape of the meniscus. Hence, for $|x-x_{c}^{\ast }|\ll \unicode[STIX]{x1D716}_{1}$ , $\unicode[STIX]{x1D6F7}$ is given by (3.29), $d_{0}^{\pm ,2\pm }$ satisfy (3.25) and $d_{1}^{\pm ,2\pm }$ satisfy (3.30).

The two regimes ( $H\gg 1$ and $H_{x}\ll \unicode[STIX]{x1D716}_{1}^{-1}\unicode[STIX]{x1D716}_{2}^{-1}$ ) overlap for $\unicode[STIX]{x1D716}_{1}\unicode[STIX]{x1D716}_{2}^{1/2}\ll |x-x_{c}^{\ast }|\ll \unicode[STIX]{x1D716}_{1}$ . By applying a similar analysis to region V, we deduce that constant solutions for $d_{0}^{\pm ,2\pm }$ are applicable throughout the water region. Thus, our previous assertion that $d_{0}^{\pm ,2\pm }\equiv 1+10^{\unicode[STIX]{x1D70F}}(\bar{c}^{1,2}-1)$ is justified. To show that $d_{1}^{\pm ,2\pm }\equiv 0$ in the thin film, we use the fact that, in regions I and III, $d^{\pm ,2\pm }\sim 1+10^{\unicode[STIX]{x1D70F}}(c^{1,2}-1)+O(\unicode[STIX]{x1D6FF})$ , since the solution for $d^{\pm ,2\pm }$ is found using an expansion in powers of $\unicode[STIX]{x1D6FF}$ . Since $d_{1}^{\pm ,2\pm }$ appear as the $O(\unicode[STIX]{x1D6FF}^{1/2}\unicode[STIX]{x1D709}^{-1/2})\gg \unicode[STIX]{x1D6FF}$ terms in the expansions of $d^{\pm ,2\pm }$ in the thin film, we must match $d_{1}^{\pm ,2\pm }$ with zero in regions IV and V. Hence, the constant solutions $d_{1}^{\pm ,2\pm }\equiv 0$ are applicable. Finally, as (3.29) is asymptotically equivalent to (3.19) as $H\rightarrow \infty$ , equation (3.29) is valid near the clay wall throughout regions IV and V.

It is worth noting that a streaming potential does not develop in the leading-order solution since we have assumed that $Pe^{\pm ,2\pm }=O(1)$ . Therefore, the terms relating to ionic advection are $O(\unicode[STIX]{x1D709}\unicode[STIX]{x1D716}_{2},\unicode[STIX]{x1D709}\unicode[STIX]{x1D6FF})$ in equations (3.21)–(3.23), and consequently affect the potential at $O(\unicode[STIX]{x1D716}_{2}^{2},\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D6FF})\ll O(\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D709}^{-1})^{2}$ ; this result is observed by Schnitzer, Frankel & Yariv (Reference Schnitzer, Frankel and Yariv2012).

Figure 4.

Schematic diagram of region IV. For $|x-x_{c}^{\ast }|\gg \unicode[STIX]{x1D716}_{1}\unicode[STIX]{x1D716}_{2}^{1/2}$ the bulk water solution is valid, as $H\gg 1$ . For $|x-x_{c}^{\ast }|\ll \unicode[STIX]{x1D716}_{1}$ the thin-film solution is valid, as $H_{x}\gg \unicode[STIX]{x1D716}_{1}^{-1}\unicode[STIX]{x1D716}_{2}^{-1}$ .

4.1.1 Shape of the oil–water interface in the transition region

In the transition region, we expect the electrostatic, capillary and viscous forces to balance. In order to balance the terms in (2.52), (2.53), (2.60) and (2.61), we define the rescalings

(4.1a-c ) $$\begin{eqnarray}\displaystyle x-(x_{c}^{\ast }+L_{o}^{\ast })=\unicode[STIX]{x1D716}_{2}^{1/2}\bar{X},\quad z=\unicode[STIX]{x1D716}_{2}Z,\quad h=\unicode[STIX]{x1D716}_{2}H, & & \displaystyle\end{eqnarray}$$

(4.2a-e ) $$\begin{eqnarray}\displaystyle u_{w}=\unicode[STIX]{x1D716}_{2}^{3/2}\bar{U}_{w},\quad u_{o}=\unicode[STIX]{x1D716}_{2}^{3/2}\bar{U}_{o},\quad w_{w}=\unicode[STIX]{x1D716}_{2}^{2}\bar{W}_{w},\quad w_{o}=\unicode[STIX]{x1D716}_{2}^{2}\bar{W}_{o},\quad t=\unicode[STIX]{x1D716}_{2}^{-1}\bar{T}, & & \displaystyle\end{eqnarray}$$

where $L_{o}^{\ast }=L_{o}/L$ is the dimensionless length of the oil slug. At leading order in the new variables, equations (2.52), (2.53), (2.60) and (2.61) read

(4.3a,b ) $$\begin{eqnarray}\displaystyle \bar{U}_{w_{ZZ}}={p_{w}}_{\bar{X}}-\frac{1}{P^{\ast }}\unicode[STIX]{x1D6F7}_{\bar{X}}\unicode[STIX]{x1D6F7}_{ZZ},\quad \bar{U}_{w_{\bar{X}}}+\bar{W}_{w_{Z}}=0, & & \displaystyle\end{eqnarray}$$

(4.4a,b ) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D707}^{\ast }\bar{U}_{o_{ZZ}}={p_{o}}_{\bar{X}},\quad \bar{U}_{o_{\bar{X}}}+\bar{W}_{o_{Z}}=0, & & \displaystyle\end{eqnarray}$$

(4.5a,b ) $$\begin{eqnarray}\displaystyle {p_{w}}_{Z}=\frac{1}{P^{\ast }}\unicode[STIX]{x1D6F7}_{Z}\unicode[STIX]{x1D6F7}_{ZZ},\quad {p_{o}}_{Z}=0, & & \displaystyle\end{eqnarray}$$

subject to

(4.6a,b ) $$\begin{eqnarray}\displaystyle \bar{U}_{w}=\bar{W}_{w}=0,\quad Z=0, & & \displaystyle\end{eqnarray}$$

(4.7a-d ) $$\begin{eqnarray}\displaystyle \bar{U}_{w}=\bar{U}_{o},\quad \bar{W}_{w}=\bar{W}_{o}=\bar{U}H_{\bar{X}}+H_{\bar{T}},\quad \bar{U}_{w_{Z}}=\unicode[STIX]{x1D707}^{\ast }\bar{U}_{o_{Z}},\quad Z=H, & & \displaystyle\end{eqnarray}$$

(4.8a,b ) $$\begin{eqnarray}\displaystyle p_{w}=p_{o}+\frac{1}{2P^{\ast }}\unicode[STIX]{x1D6F7}_{Z}^{2}|_{Z=H^{-}}-\frac{\unicode[STIX]{x1D6E4}}{2P^{\ast }}H_{\bar{X}\bar{X}},\quad Z=H, & & \displaystyle\end{eqnarray}$$

(4.9a-c ) $$\begin{eqnarray}\displaystyle {U_{o}}_{Z}\rightarrow 0,\quad {W_{o}}_{Z}\rightarrow 0,\quad Z\rightarrow \infty , & & \displaystyle\end{eqnarray}$$

where we have substituted $\unicode[STIX]{x1D6F7}_{Z}=0$ in the oil phase in (4.8). Note that (4.9) is the constraint that the velocity is bounded in the bulk of the oil.

Integrating (4.5b ) in the oil phase yields $p_{o}=p_{o}(\bar{X},\bar{T})$ . By substituting this into (4.4a ) and applying (4.9), we find that $\bar{U}_{o_{Z}}=0$ and $p_{o_{\bar{X}}}=0$ . Hence, $\bar{U}_{o}=\bar{U}_{o}(\bar{X},\bar{T})$ and $p_{o}=p_{o}(\bar{T})$ . By integrating (4.5a ) and applying (4.8), we find that the pressure in the water layer is given by

(4.10) $$\begin{eqnarray}\displaystyle p_{w}=p_{o}+\frac{1}{2P^{\ast }}\unicode[STIX]{x1D6F7}_{Z}^{2}-\frac{\unicode[STIX]{x1D6E4}}{2P^{\ast }}H_{\bar{X}\bar{X}}. & & \displaystyle\end{eqnarray}$$

By substituting this into (4.3a ), and applying (4.6) and (4.7), we find that the horizontal velocity in the transition region is given by

(4.11) $$\begin{eqnarray}\displaystyle \bar{U}_{w}=\frac{1}{4P^{\ast }}(2HZ-Z^{2})\left(\unicode[STIX]{x1D6E4}H_{\bar{X}\bar{X}}+\frac{{\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}\cosh H+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2}}{\sinh ^{2}H}\right)_{\bar{X}}, & & \displaystyle\end{eqnarray}$$

where we have substituted $\unicode[STIX]{x1D6F7}$ using (3.34). By integrating (4.3b ) from $Z=0$ to $Z=H$ , and applying (4.6) and (4.7), we obtain the following thin-film equation describing the thickness of the water layer, $H(\bar{X},\bar{T})$ , in the transition region:

(4.12) $$\begin{eqnarray}\displaystyle H_{\bar{T}}+\left[\frac{H^{3}}{6P^{\ast }}\left(\unicode[STIX]{x1D6E4}H_{\bar{X}\bar{X}}+\frac{{\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}\cosh H+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2}}{\sinh ^{2}H}\right)_{\bar{X}}\right]_{\bar{X}}=0. & & \displaystyle\end{eqnarray}$$

In the limit $\bar{X}\rightarrow -\infty$ we impose the condition that the solution approaches a constant value, $H^{\ast }$ , to be determined. As $\bar{X}\rightarrow \infty$ , we require the solution to match with the meniscus region. We assume that the radius of the meniscus is approximately the radius of the capillary tube and that the leading-order correction, $r_{1}^{+}$ , is of the order of the thickness of the wetting layer, i.e. $r^{+}=1-\unicode[STIX]{x1D716}_{2}r_{1}^{+}+O(\unicode[STIX]{x1D716}_{2}^{2})$ (the $+$ superscripts denote the front meniscus). Hence, the boundary conditions we impose on $H$ are $H\rightarrow H^{\ast }$ as $\bar{X}\rightarrow -\infty$ and $H_{\bar{X}\bar{X}}\rightarrow 1$ as $\bar{X}\rightarrow \infty$ .

We assume a quasi-static travelling wave solution, $H(\bar{X},\bar{T})=H(Y)$ , where $Y=\bar{X}-\bar{X}_{R}(\bar{T})$ , and $\bar{X}_{R}$ corresponds with the origin in the reference frame moving with the oil slug, which can be chosen arbitrarily. In the travelling wave coordinate, equation (4.12) becomes

(4.13) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }VH_{Y}=\left[\frac{H^{3}}{6}\left(\unicode[STIX]{x1D6E4}H_{YY}+\frac{{\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}\cosh H+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2}}{\sinh ^{2}H}\right)_{Y}\right]_{Y}, & & \displaystyle\end{eqnarray}$$

subject to

(4.14) $$\begin{eqnarray}\displaystyle & \displaystyle H\rightarrow H^{\ast },\quad \text{as }Y\rightarrow -\infty , & \displaystyle\end{eqnarray}$$

(4.15) $$\begin{eqnarray}\displaystyle & \displaystyle H_{YY}\rightarrow 1,\quad \text{as }Y\rightarrow \infty , & \displaystyle\end{eqnarray}$$

where $V=\unicode[STIX]{x1D716}_{2}^{3/2}\,\text{d}\bar{X}_{R}/\text{d}\bar{T}$ is the velocity at which the oil slug moves, scaled with the macroscopic variables. (In appendix A, we show that (4.14) and (4.15) impose the required number of constraints on (4.16) to be able to find both $H(Y)$ and $H^{\ast }$ .) We integrate (4.13) once to arrive at

(4.16) $$\begin{eqnarray}\displaystyle 6\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }V\frac{(H-H^{\ast })}{H^{3}}=\unicode[STIX]{x1D6E4}H_{YYY}+\left(\frac{{\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}\cosh H+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2}}{\sinh ^{2}H}\right)_{Y}. & & \displaystyle\end{eqnarray}$$

We would like to be able to determine the film thickness, $H^{\ast }$ , for a given applied pressure, $P^{\ast }$ . However, at this stage, the velocity, $V$ , is still unknown. To better understand (4.16), we will, for the moment, treat $P^{\ast }V$ as an independent parameter. In § 4.3 we will find $V$ by solving the macroscopic flow problem in combination with (4.16). We note that we can rescale the transition region equation, (4.16), so that the viscous and capillary forces balance and we recover a Bretherton-type equation (Bretherton Reference Bretherton1961), which includes electrostatic effects. The details of this rescaling are given in appendix B.

To solve (4.16) subject to (4.14) and (4.15) for a given $P^{\ast }V$ , we use a shooting method. We treat (4.16) as an initial value problem with $H^{\ast }$ as the initial parameter, and vary $H^{\ast }$ until $H(Y)$ satisfies (4.15) in the limit $Y\rightarrow \infty$ . This gives us the correct value of the film thickness, $H^{\ast }$ , and the shape of the transition region, $H(Y)$ .

In figure 7(b) we plot $h(x)=\unicode[STIX]{x1D716}_{2}H(Y)$ in the high-salinity ( $\unicode[STIX]{x1D70F}=0$ ) and low-salinity ( $\unicode[STIX]{x1D70F}=-6$ ) regimes, for $P^{\ast }V=9\times 10^{-6}$ , where $P^{\ast }$ is determined by the values in table 2, and we have chosen $V=1$ . Values for the other parameters used are given in table 1. In the numerical computations, we use a fourth-order Runge–Kutta method to solve (4.16). We observe that the film is thicker in the low-salinity regime than in the high-salinity regime. This is consistent with the MIE proposal that reducing the salinity of the brine increases the film thickness, thus aiding the release of oil.

Figure 7.

Plots of the rear (a) and front (b) transition region profiles, $h(x)$ , in the low-salinity ( $\unicode[STIX]{x1D70F}=-6$ ) and high-salinity ( $\unicode[STIX]{x1D70F}=0$ ) regimes, for $P^{\ast }V=9\times 10^{-6}$ .

4.1.2 Correction to the radius of curvature

To determine the leading-order correction to the radius of curvature in the meniscus region, $r_{1}^{+}$ , we consider the behaviour of $H(Y)$ as $Y$ approaches the meniscus region. Since $H_{YY}\sim 1$ as $Y\rightarrow \infty$ , at leading order

(4.17) $$\begin{eqnarray}\displaystyle H\sim {\textstyle \frac{1}{2}}(Y+\unicode[STIX]{x1D6FC}^{+})^{2}+\unicode[STIX]{x1D6FD}^{+}, & & \displaystyle\end{eqnarray}$$

as $Y\rightarrow \infty$ for some constants $\unicode[STIX]{x1D6FC}^{+}$ and $\unicode[STIX]{x1D6FD}^{+}$ . Using the leading-order translational invariance of the problem, $\unicode[STIX]{x1D6FC}^{+}$ can be chosen arbitrarily. To determine $\unicode[STIX]{x1D6FD}^{+}$ , we match (4.17) with the film thickness in the meniscus region. In the meniscus region the oil–water interface is circular, which, in terms of the macroscopic variables $x$ and $h$ , gives

(4.18) $$\begin{eqnarray}\displaystyle (1-h)^{2}+(x-x^{\ast })^{2}=(r^{+})^{2}\sim (1-\unicode[STIX]{x1D716}_{2}r_{1}^{+})^{2}, & & \displaystyle\end{eqnarray}$$

where $(x^{\ast },1)$ is the centre of the circle with the same tangent and curvature as the extremal point of the oil–water interface, as shown in figure 8. As $x\rightarrow x^{\ast }$ (i.e. as $x$ approaches the transition region), expanding in $(x-x^{\ast })$ and $\unicode[STIX]{x1D716}_{2}$ yields

(4.19) $$\begin{eqnarray}\displaystyle h\sim {\textstyle \frac{1}{2}}(x-x^{\ast })^{2}+\unicode[STIX]{x1D716}_{2}r_{1}^{+}+O((x-x^{\ast })^{4},\unicode[STIX]{x1D716}_{2}(x-x^{\ast })^{2}). & & \displaystyle\end{eqnarray}$$

Matching the inner solution (4.17) with the outer solution (4.19) gives $\unicode[STIX]{x1D6FD}^{+}=r_{1}^{+}$ . Rearranging (4.17) then gives the first correction to the radius of the meniscus as

(4.20) $$\begin{eqnarray}\displaystyle r_{1}^{+}=\lim _{Y\rightarrow \infty }H-{\textstyle \frac{1}{2}}H_{Y}^{2}. & & \displaystyle\end{eqnarray}$$

By substituting the solution for $H(Y)$ found in § 4.1.1 into (4.20), we can therefore compute the correction to the curvature in the front meniscus region.

Figure 8.

Schematic diagram of the circular solution the oil–water interface approaches.

4.2.1 Shape of the oil–water interface in the transition region

We study the rear end of the oil slug advancing with speed $V$ by considering the equivalent problem of the front end receding with speed $V$ . Thus, to find the shape of the oil–water interface at the rear end of the oil slug, we make the substitution $V\mapsto -V$ and solve (4.16) subject to (4.14) and (4.15) using a shooting method. However, as we found $H^{\ast }$ by considering the shape of the front transition region, we can no longer use this as a shooting parameter. To determine the correct parameter to vary, we examine the behaviour of $H(Y)$ as $Y\rightarrow -\infty$ (towards the thin-film region). In this limit, we assume that $H\sim H^{\ast }+\bar{H}$ , where $\bar{H}$ is small. By substituting this into (4.16), applying $\bar{H}\rightarrow 0$ as $Y\rightarrow -\infty$ , and using the translational invariance of the problem to eliminate one constant of integration, we find that

(4.21) $$\begin{eqnarray}\displaystyle H\sim H^{\ast }+\text{e}^{k_{1}Y}+A\text{e}^{k_{2}Y}. & & \displaystyle\end{eqnarray}$$

The constants $k_{1}$ and $k_{2}$ are roots of the characteristic equation associated with the linear ordinary differential equation (ODE) satisfied by $\bar{H}$ (see appendix A for details), and $A$ is an unknown constant. We determine the shape of the rear transition region, $H(Y)$ , by solving (4.16) (with the substitution $V\mapsto -V$ ) subject to (4.15) in the limit $Y\rightarrow \infty$ and (4.21) in the limit $Y\rightarrow -\infty$ , with $A$ as a shooting parameter.

In figure 7(a) we plot $h(x)=\unicode[STIX]{x1D716}_{2}H(Y)$ in the high-salinity ( $\unicode[STIX]{x1D70F}=0$ ) and low-salinity ( $\unicode[STIX]{x1D70F}=-6$ ) regimes, with $P^{\ast }V=9\times 10^{-6}$ . We observe that the profile of the rear of the oil slug is non-monotonic in the high-salinity regime. This is in contrast to the profile of the front of the oil slug, shown in figure 7(b), which is monotonic in both regimes. This is consistent with Bretherton (Reference Bretherton1961), who studied the motion of a bubble through a capillary and found that the front meniscus is monotonic, whereas the rear meniscus exhibits sinusoidal behaviour in the absence of electrostatic forces.

4.2.2 Correction to the radius of curvature

At the rear of the oil slug, we find that, similarly to the front of the oil slug, the leading-order correction to the radius of curvature of the meniscus is given by

(4.22) $$\begin{eqnarray}\displaystyle r_{1}^{-}=\lim _{Y\rightarrow \infty }H-{\textstyle \frac{1}{2}}H_{Y}^{2}, & & \displaystyle\end{eqnarray}$$

where the radius of curvature is $r^{-}\sim 1-\unicode[STIX]{x1D716}_{2}r_{1}^{-}+O(\unicode[STIX]{x1D716}_{2}^{2})$ , and $H(Y)$ is the solution found in § 4.2.1.

4.3.1 Regions I and III

From (4.23b ), we see that $p_{w}=p_{w}(x)$ . By integrating (4.23a ) twice and applying (4.29) and (4.31), we find that, in regions I and III,

(4.32) $$\begin{eqnarray}\displaystyle u_{w}={p_{w}}_{x}({\textstyle \frac{1}{2}}z^{2}-z). & & \displaystyle\end{eqnarray}$$

By integrating (4.23c ) from $z=0$ to $z=1$ and applying (4.29) and (4.31), we find that

(4.33) $$\begin{eqnarray}\displaystyle \frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x}\left[\int _{0}^{1}u_{w}\,\text{d}z\right]=0. & & \displaystyle\end{eqnarray}$$

It follows that ${p_{w}}_{xx}=0$ , i.e. the pressure gradient is constant. Integrating (4.33) and substituting (4.32) yields

(4.34) $$\begin{eqnarray}\displaystyle Q_{w}:=\int _{0}^{1}u_{w}\,\text{d}z=-\frac{1}{3}{p_{w}}_{x}, & & \displaystyle\end{eqnarray}$$

where $Q_{w}$ is the volumetric flow rate of water in the lower half of the capillary tube. By applying (4.25)–(4.28) we find that

(4.35) $$\begin{eqnarray}\displaystyle & \displaystyle Q_{w}=\frac{1}{3x_{c}^{\ast }}(1-P_{1}),\quad \text{in region I,} & \displaystyle\end{eqnarray}$$

(4.36) $$\begin{eqnarray}\displaystyle & \displaystyle Q_{w}=\frac{1}{3(1-x_{c}^{\ast }-L_{o}^{\ast })}P_{2},\quad \text{in region III.} & \displaystyle\end{eqnarray}$$

4.3.2 Region II

From (4.24b ), we see that $p_{o}=p_{o}(x)$ . Integrating (4.23a ) and (4.24a ) subject to (4.29)–(4.31) yields (to $O(\unicode[STIX]{x1D716}_{2})$ )

(4.37) $$\begin{eqnarray}\displaystyle & \displaystyle u_{w}={p_{w}}_{x}({\textstyle \frac{1}{2}}z^{2}-z),\quad 0<z<\unicode[STIX]{x1D716}_{2}H^{\ast }, & \displaystyle\end{eqnarray}$$

(4.38) $$\begin{eqnarray}\displaystyle & \displaystyle u_{o}=\frac{1}{\unicode[STIX]{x1D707}^{\ast }}{p_{o}}_{x}\left(\frac{1}{2}z^{2}-z+\unicode[STIX]{x1D716}_{2}H^{\ast }(1-\unicode[STIX]{x1D707}^{\ast })\right),\quad \unicode[STIX]{x1D716}_{2}H^{\ast }<z<1. & \displaystyle\end{eqnarray}$$

By integrating (4.23c ) from $z=0$ to $z=\unicode[STIX]{x1D716}_{2}H^{\ast }$ , (4.24c ) from $z=\unicode[STIX]{x1D716}_{2}H^{\ast }$ to $z=1$ , and applying (4.29)–(4.31), we find that

(4.39) $$\begin{eqnarray}\displaystyle \frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x}\left[\int _{0}^{\unicode[STIX]{x1D716}_{2}H^{\ast }}u_{w}\,\text{d}z+\int _{\unicode[STIX]{x1D716}_{2}H^{\ast }}^{1}u_{o}\,\text{d}z\right]=0,\quad \text{in region II}. & & \displaystyle\end{eqnarray}$$

By substituting (4.37), we find that the water flux is $O(\unicode[STIX]{x1D716}_{2}^{2})$ , since the water layer is thin and has low velocity near the capillary surface, and can thus be neglected. Hence, by substituting (4.38), it follows that ${p_{o}}_{xx}=0$ . Integrating (4.39) yields

(4.40) $$\begin{eqnarray}\displaystyle Q_{o}:=\int _{\unicode[STIX]{x1D716}_{2}H}^{1}u_{o}\,\text{d}z=-\frac{1}{\unicode[STIX]{x1D707}^{\ast }}{p_{o}}_{x}\left(\frac{1}{3}-\unicode[STIX]{x1D716}_{2}(1-\unicode[STIX]{x1D707}^{\ast })H^{\ast }\right),\quad \text{in region II}, & & \displaystyle\end{eqnarray}$$

where $Q_{o}$ is the volumetric flow rate of oil in the lower half of the capillary tube. By applying (4.26) and (4.27), we find that

(4.41) $$\begin{eqnarray}\displaystyle Q_{o}=\frac{1}{\unicode[STIX]{x1D707}^{\ast }L_{o}^{\ast }}\left\{\frac{P_{1}-P_{2}}{3}+\unicode[STIX]{x1D716}_{2}\left[\frac{\unicode[STIX]{x1D6E4}}{6P^{\ast }}(r_{1}^{-}-r_{1}^{+})-(P_{1}-P_{2})(1-\unicode[STIX]{x1D707}^{\ast })H^{\ast }\right]\right\},\quad \text{in region II.} & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

To better understand how the water film affects the flow of oil, we derive the effective boundary condition satisfied by the oil on the capillary surface. By evaluating $u_{o}$ and ${u_{o}}_{z}$ on $z=0$ using (4.38), and eliminating ${p_{o}}_{x}$ , we find that

(4.42) $$\begin{eqnarray}\displaystyle u_{o}=\unicode[STIX]{x1D716}_{2}H^{\ast }(\unicode[STIX]{x1D707}^{\ast }-1){u_{o}}_{z},\quad \text{on }z=0. & & \displaystyle\end{eqnarray}$$

We observe that, since the oil is more viscous than the water, the presence of the thin water film below the oil slug results in partial slip along the capillary surface, with the slip length being proportional to the width of the water film.

4.3.3 Macroscopic flux

Since we assume that the slug fills almost the whole width of the channel ( $\unicode[STIX]{x1D716}_{2}H^{\ast }$ is small), conservation of volume imposes that $Q:=Q_{o}\sim Q_{w}$ . Hence, by solving the three equations, (4.35), (4.36) and (4.41), for the three unknowns, $P_{1}$ , $P_{2}$ and $Q$ , we find that

(4.43) $$\begin{eqnarray}\displaystyle Q\sim \frac{1}{1+(\unicode[STIX]{x1D707}^{\ast }-1)L_{o}^{\ast }}\left\{\frac{1}{3}+\unicode[STIX]{x1D716}_{2}\left[\frac{\unicode[STIX]{x1D6E4}}{6P^{\ast }}(r_{1}^{-}-r_{1}^{+})+\frac{\unicode[STIX]{x1D707}^{\ast }L_{o}^{\ast }}{1+(\unicode[STIX]{x1D707}^{\ast }-1)L_{o}^{\ast }}(\unicode[STIX]{x1D707}^{\ast }-1)H^{\ast }\right]\!\right\}. & & \displaystyle\end{eqnarray}$$

We see from (4.43) that the volumetric flow rate consists of a leading-order term, plus a capillary contribution resulting from the difference between the droplet shape at the front and rear of the slug, and a slip contribution resulting from the finite film thickness. In order to make use of this expression, we require the solutions for $r_{1}^{\pm }$ and $H^{\ast }$ found in §§ 4.1 and 4.2 by solving (4.16) subject to (4.14) and (4.15). However, since the volumetric flow rate is the velocity of the oil slug multiplied by the channel width, $Q=V$ , the solution for $H$ given by solving (4.16) is dependent on $Q$ . Thus, we must solve (4.16) and (4.43) simultaneously to determine the velocity.

In figure 10(a), we plot the film thickness, $h^{\ast }=\unicode[STIX]{x1D716}_{2}H^{\ast }$ , against the salinity, $\unicode[STIX]{x1D70F}$ . As the salinity increases from the low-salinity regime to the high-salinity regime, the film thickness decreases by 39 %. This is consistent with the MIE model proposed by Lager et al. (Reference Lager, Webb, Black, Singleton and Sorbie2008). However, the percentage decrease is larger than found experimentally by Lee et al. (Reference Lee, Webb, Collins, Lager, Clarke, O’Sullivan, Routh and Wang2010). This is probably due to the fact that structural and van der Waals forces have been neglected in this model. In the high-salinity limit, structural forces dominate the electrostatic forces, which prevents the film thickness from reducing further. Hirasaki (Reference Hirasaki1991) discusses this in more detail. Furthermore, as the thickness of the film becomes comparable with the size of the ions in the model, the continuum assumption breaks down and the model presented here is no longer applicable.

In figure 10(b) we plot the corrections to the menisci radii, $r_{1}^{\pm }$ , and the magnitude of the difference between them, $|r_{1}^{+}-r_{1}^{-}|$ . The front and rear corrections are similar, so that the difference between them is an order of magnitude smaller than the corrections themselves. The difference is approximately 8.4 times larger in the high-salinity regime than in the low-salinity regime.

Figure 10.

(a) Plot of the film thickness, $h^{\ast }$ , as a function of the salinity, $\unicode[STIX]{x1D70F}$ . (b) Plots of the corrections to the menisci radii, $r_{1}^{\pm }$ , and the difference between the front and rear menisci, $|r_{1}^{+}-r_{1}^{-}|$ , as functions of the salinity, $\unicode[STIX]{x1D70F}$ . The low-salinity and high-salinity regimes correspond to $\unicode[STIX]{x1D70F}\approx -6$ , and $\unicode[STIX]{x1D70F}\approx 0$ , respectively.

In figure 11(a), we plot the velocity of the oil slug, $V$ , against the salinity, $\unicode[STIX]{x1D70F}$ . We choose $L_{o}^{\ast }=0.5$ for this calculation and for subsequent computations. We observe that the velocity decreases as the salinity increases, with the velocity being approximately 3.1 times larger in the low-salinity regime than in the high-salinity regime.

In figure 11(b), we plot the components of the velocity corresponding to (i) the difference in the menisci radii, $V_{1}^{(1)}=\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D6E4}(r_{1}^{-}-r_{1}^{+})/6P^{\ast }$ , and (ii) the slip resulting from the thin water layer, $V_{1}^{(2)}=\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D707}^{\ast }L_{o}^{\ast }(\unicode[STIX]{x1D707}^{\ast }-1)H^{\ast }/(1+(\unicode[STIX]{x1D707}^{\ast }-1)L_{o}^{\ast })$ . We observe that these components are larger in magnitude than the total velocity, $V$ . Since $P^{\ast }\ll \unicode[STIX]{x1D716}_{2}$ for the values of the parameters used (given in table 2), we find that the primary balance in (4.43) is between $V_{1}^{(1)}$ (which depends on $V$ through (4.16)) and $1/3+V_{1}^{(2)}$ . This is in contrast with the $P^{\ast }\gg \unicode[STIX]{x1D716}_{2}$ regime, in which $V_{1}^{(1)}$ and $V_{1}^{(2)}$ are perturbations to the leading-order velocity. Although $P^{\ast }\ll \unicode[STIX]{x1D716}_{2}$ for the values used in this paper, it is worth noting that, due to the diversity of reservoir conditions, there will probably be applications in which $P^{\ast }\gg \unicode[STIX]{x1D716}_{2}$ . Since we are able to make further analytical progress in this regime, we will briefly discuss the case $P^{\ast }\gg \unicode[STIX]{x1D716}_{2}$ in § 4.4.

Figure 11.

Plots of the velocity of the oil slug, given by (4.43), as functions of the salinity, $\unicode[STIX]{x1D70F}$ , where $\unicode[STIX]{x1D70F}\approx -6$ corresponds to the low-salinity regime, and $\unicode[STIX]{x1D70F}\approx 0$ corresponds to the high-salinity regime: (a) the total velocity, and (b) the two corrections to the velocity. The latter are $V_{1}^{(1)}=(\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D6E4}/6P^{\ast })(r_{1}^{-}-r_{1}^{+})$ and $V_{1}^{(2)}=(\unicode[STIX]{x1D716}_{2}\unicode[STIX]{x1D707}^{\ast }L_{o}^{\ast }/(1+(\unicode[STIX]{x1D707}^{\ast }-1)L_{o}^{\ast }))(\unicode[STIX]{x1D707}^{\ast }-1)H^{\ast }$ .

4.5.1 Film thickness in the limit $P^{\ast }\ll \unicode[STIX]{x1D716}_{2}^{3/2}$

To determine the film thickness in the limit $P^{\ast }\ll \unicode[STIX]{x1D716}_{2}^{3/2}$ , we return to (4.16) and integrate once to find that

(4.45) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6E4}H_{YY}+\frac{{\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}\cosh H+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2}}{\sinh ^{2}H}+6\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }V\int _{Y}^{\infty }\frac{H-H^{\ast }}{H^{3}}\,\text{d}Y^{\prime }=\unicode[STIX]{x1D6E4}. & & \displaystyle\end{eqnarray}$$

We expand $H\sim H_{0}+\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }H_{1}+O(\unicode[STIX]{x1D716}_{2}^{-3}{P^{\ast }}^{2})$ and $H^{\ast }\sim H_{0}^{\ast }+\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }H_{1}^{\ast }+O(\unicode[STIX]{x1D716}_{2}^{-3}{P^{\ast }}^{2})$ and consider (4.45) in the limit $Y\rightarrow -\infty$ . We find that, at leading order, after some simplification, the film thickness is given by

(4.46) $$\begin{eqnarray}\displaystyle H_{0}^{\ast }=\cosh ^{-1}\!\left(\frac{s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}+\sqrt{({\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+\unicode[STIX]{x1D6E4})({s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2}+\unicode[STIX]{x1D6E4})}}{\unicode[STIX]{x1D6E4}}\right). & & \displaystyle\end{eqnarray}$$

In figure 12(a) we plot the film thickness against the salinity, for $-6<\unicode[STIX]{x1D70F}<0$ . We compare the numerical solution found in § 4.3.3, using $P^{\ast }=9\times 10^{-6}$ , with the asymptotic solution in the limit $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }\rightarrow 0$ , equation (4.46). There is good agreement between the two solutions. In figure 12(b) we show that the film thickness, $h^{\ast }=\unicode[STIX]{x1D716}_{2}H^{\ast }$ , found numerically in § 4.3.3 varies with $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }$ in the low-salinity regime, for $10^{-2}<\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }<10^{2}$ . We compare the result with the film thickness found in the asymptotic limit $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }\rightarrow 0$ . The asymptotic solution agrees with the numerical solution to within 4.4 %, for $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }<10$ .

Figure 12.

(a) Film thickness as a function of the salinity, comparing the numerical solution for $P^{\ast }=9\times 10^{-6}$ with the asymptotic solution in the limit $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }\rightarrow 0$ . (b) Film thickness as a function of the dimensionless pressure, $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }$ , in the low-salinity regime, comparing the numerical solution with the analytic solution in the limit $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }\rightarrow 0$ .

4.5.2 Corrections to the menisci radii in the limit $P^{\ast }\ll \unicode[STIX]{x1D716}_{2}^{3/2}$

To determine $r_{1}^{\pm }$ in the quasi-static limit, we integrate (4.45) to get

(4.47) $$\begin{eqnarray}\displaystyle & & \displaystyle \frac{\unicode[STIX]{x1D6E4}}{2}H_{Y}^{2}+\frac{2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}+({\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2})\cosh H^{\ast }}{\sinh H^{\ast }}-\frac{2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}+({\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2})\cosh H}{\sinh H}\nonumber\\ \displaystyle & & \displaystyle \qquad +\,6\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }V\int _{-\infty }^{Y}H_{Y}\int _{Y^{\prime }}^{\infty }\frac{H-H^{\ast }}{H^{3}}\,\text{d}Y^{\prime \prime }\,\text{d}Y^{\prime }=\unicode[STIX]{x1D6E4}(H-H^{\ast }).\end{eqnarray}$$

We expand $H\sim H_{0}+\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }H_{1}+O(\unicode[STIX]{x1D716}_{2}^{-3}{P^{\ast }}^{2})$ and consider the leading-order terms in the limit $Y\rightarrow \infty$ to get

(4.48) $$\begin{eqnarray}\displaystyle H_{0}\sim \frac{1}{2}(Y+\unicode[STIX]{x1D6FC}^{+})^{2}+\left\{H_{0}^{\ast }+\frac{({\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2})}{\unicode[STIX]{x1D6E4}}(\coth H_{0}^{\ast }-1)+\frac{2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}}{\unicode[STIX]{x1D6E4}\sinh H_{0}^{\ast }}\right\}. & & \displaystyle\end{eqnarray}$$

By matching terms with (4.17), and using the fact that $\unicode[STIX]{x1D6FD}^{+}=r_{1}^{+}$ , we find that the bracketed term in (4.48) must be the first-order (in $\unicode[STIX]{x1D716}_{2}$ ) correction to the radius of the meniscus,

(4.49) $$\begin{eqnarray}\displaystyle r_{1}^{+}\sim H_{0}^{\ast }+\frac{({\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2})}{\unicode[STIX]{x1D6E4}}(\coth H_{0}^{\ast }-1)+\frac{2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}}{\unicode[STIX]{x1D6E4}\sinh H_{0}^{\ast }}+O(\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }). & & \displaystyle\end{eqnarray}$$

Since the oil slug is moving, we expect the front and rear menisci to have different radii, and hence, by the Laplace–Young equation (2.61), the pressure drop across the oil will be dependent on the difference between these two radii. To find the correction to the meniscus radius at the rear end of the oil slug, we make the substitution $V\mapsto -V$ . However, as the leading-order correction to the radius of curvature (4.49) is independent of $V$ , the first-order correction is required in order to find the leading-order difference between the radii of the two menisci.

The $O(\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast })$ terms in equation (4.47) read

(4.50) $$\begin{eqnarray}\displaystyle & & \displaystyle \unicode[STIX]{x1D6E4}{H_{0}}_{Y}{H_{1}}_{Y}+H_{1}\frac{2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}\cosh H_{0}+{\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2}}{\sinh ^{2}H_{0}}\nonumber\\ \displaystyle & & \displaystyle \qquad +\,6V\int _{-\infty }^{Y}{H_{0}}_{Y}\int _{Y^{\prime }}^{\infty }\frac{H_{0}-H_{0}^{\ast }}{H_{0}^{3}}\,\text{d}Y^{\prime \prime }\,\text{d}Y^{\prime }=\unicode[STIX]{x1D6E4}H_{1},\end{eqnarray}$$

where we have made use of (4.46) in eliminating some of the terms. We find that

(4.51) $$\begin{eqnarray}\displaystyle H_{1}\sim \frac{6V}{\unicode[STIX]{x1D6E4}}\int _{-\infty }^{\infty }\frac{(H_{0}-H_{0}^{\ast })^{2}}{H_{0}^{3}}\,\text{d}Y^{\prime }, & & \displaystyle\end{eqnarray}$$

as $Y\rightarrow \infty$ , which provides the $O(\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast })$ term in the expression for $r_{1}^{+}$ . By making the substitution $V\mapsto -V$ , we find that the $O(\unicode[STIX]{x1D716}_{2})$ corrections to the front $(r^{+})$ and back $(r^{-})$ radii are

(4.52) $$\begin{eqnarray}\displaystyle r_{1}^{\pm }\sim H_{0}^{\ast }+\frac{({\unicode[STIX]{x1D70E}_{c}}_{0}^{2}+{s^{\ast }}^{2}{\unicode[STIX]{x1D70E}_{o}}_{0}^{2})}{\unicode[STIX]{x1D6E4}}(\coth H_{0}^{\ast }-1)+\frac{2s^{\ast }{\unicode[STIX]{x1D70E}_{c}}_{0}{\unicode[STIX]{x1D70E}_{o}}_{0}}{\unicode[STIX]{x1D6E4}\sinh H_{0}^{\ast }}\pm \frac{6P^{\ast }V}{\unicode[STIX]{x1D716}_{2}^{3/2}\unicode[STIX]{x1D6E4}}\int _{-\infty }^{\infty }\frac{(H_{0}-H_{0}^{\ast })^{2}}{H_{0}^{3}}\,\text{d}Y^{\prime }. & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

4.5.3 Velocity of the oil slug in the limit $P^{\ast }\ll \unicode[STIX]{x1D716}_{2}^{3/2}$

To determine the velocity of the oil slug in the limit $P^{\ast }\ll \unicode[STIX]{x1D716}_{2}^{3/2}$ , we substitute (4.52) into (4.43) and rearrange (using the fact that $Q=V$ ) to find

(4.53) $$\begin{eqnarray}\displaystyle Q=V\sim \frac{\displaystyle \frac{1+(\unicode[STIX]{x1D707}^{\ast }-1)L_{o}^{\ast }}{3}+\unicode[STIX]{x1D716}_{2}[L_{o}^{\ast }\unicode[STIX]{x1D707}^{\ast }(\unicode[STIX]{x1D707}^{\ast }-1)H_{0}^{\ast }]}{[1+(\unicode[STIX]{x1D707}^{\ast }-1)L_{o}^{\ast }]^{2}+2\unicode[STIX]{x1D716}_{2}^{-1/2}(1+(\unicode[STIX]{x1D707}^{\ast }-1)L_{o}^{\ast })\displaystyle \int _{-\infty }^{\infty }\frac{(H_{0}-H_{0}^{\ast })^{2}}{H_{0}^{3}}\,\text{d}Y^{\prime }}, & & \displaystyle\end{eqnarray}$$

where $H_{0}^{\ast }$ is given by (4.46).

It should be noted that (4.16) needs to be solved numerically to evaluate the integral term in (4.53), even in this asymptotic limit. However, this computation can be made significantly simpler. Since the integral term occurs in the first-order (in $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }$ ) correction to the meniscus radius, equation (4.52), we only require the leading-order term of the integral (as we have neglected $O(\unicode[STIX]{x1D716}_{2}^{-3}{P^{\ast }}^{2})$ terms). To determine the leading-order film profile, $H_{0}(Y)$ , we set the left-hand side of (4.16) to zero and integrate twice to obtain (4.47) (with $P^{\ast }$ set to zero), and thus we only need to solve a first-order ODE in the small $P^{\ast }$ limit. Further, since we have an expression for the film thickness, equation (4.46), we can solve (4.47) as a boundary value problem, whereas for the $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }=O(1)$ problem we need to use a shooting method, as this boundary condition is unknown.

In figure 13(a) we plot the velocity of the oil slug against the salinity, $\unicode[STIX]{x1D70F}$ , for $-6<\unicode[STIX]{x1D70F}<0$ . We compare the numerical solution found in § 4.3.3 using $P^{\ast }=9\times 10^{-6}$ , with the asymptotic solution in the limit $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }\rightarrow 0$ , equation (4.53). There is good agreement between the two solutions. The velocity is approximately three times greater in the low-salinity regime, for $P^{\ast }=0$ , which suggests that the MIE mechanism, as described by the model presented here, contributes to the low-salinity effect.

Figure 13.

(a) Velocity of the oil slug as a function of the salinity, comparing the numerical solution for $P^{\ast }=9\times 10^{-6}$ with the asymptotic solution in the limit $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }\rightarrow 0$ . (b) Velocity as a function of the dimensionless pressure, $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }$ , in the low-salinity regime, comparing the numerical solution with the analytic solution in the limit $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }\rightarrow 0$ .

In figure 13(b) we plot the velocity as a function of $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }$ in the low-salinity regime, comparing the numerical and asymptotic solutions. The solutions agree to within 2.8 % for $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }<10$ . We expect that $\unicode[STIX]{x1D716}_{2}^{-3/2}P^{\ast }<10$ for most reservoir conditions, and consequently that the quasi-static limit can usually be applied. However, there will be many applications where the applied pressure is much larger, and $P^{\ast }$ can no longer be considered small. We briefly study how larger pressures affect the shape and velocity of the oil slug.