Controlling capillary fingering using pore size gradients in disordered media

TL;DR

This study demonstrates that introducing pore size gradients in disordered porous media can effectively control capillary fingering, with experimental and computational evidence showing suppression of fingering at sufficient gradient magnitudes.

Contribution

The paper reveals how pore size gradients influence capillary fingering, providing a new method to control flow in porous media through gradient manipulation.

Findings

Pore size gradients can suppress capillary fingering.

Flow behavior depends on gradient direction and magnitude.

A non-dimensional parameter describes flow regimes.

Abstract

Capillary fingering is a displacement process that can occur when a non-wetting fluid displaces a wetting fluid from a homogeneous disordered porous medium. Here, we investigate how this process is influenced by a pore size gradient. Using microfluidic experiments and computational pore-network models, we show that the non-wetting fluid displacement behavior depends sensitively on the direction and the magnitude of the gradient. The fluid displacement depends on the competition between a pore size gradient and pore-scale disorder; indeed, a sufficiently large gradient can completely suppress capillary fingering. By analyzing capillary forces at the pore scale, we identify a non-dimensional parameter that describes the physics underlying these diverse flow behaviors. Our results thus expand the understanding of flow in complex porous media, and suggest a new way to control flow behavior via the introduction of pore size gradients.

Equations

$$V_{\text{F}}\approx {\int }_{a_{n,\text{min}}}^{a_{n,\text{max}}}{p}_n(a_n){\Phi }_{n+1}(a_n)\mathrm{d}{a}_n$$

$$V_{F}\approx\frac{1}{\sqrt{2\pi}}\int_{a_{n,\text{min}}}^{a_{n,\text{max}}}\Big{[}1+\operatorname{erf}\Big{(}\frac{\xi-g/2\sigma}{\sqrt{2}/2}\Big{)}\Big{]}e^{-2\xi^{2}}\mathrm{d}\xi$$

Full text

I Introduction

Drainage is the process by which a non-wetting fluid displaces a wetting fluid from a porous medium. This phenomenon is ubiquitous: it arises in diverse settings including groundwater contamination, oil migration, gas venting from sediments, CO2 sequestration, mercury porosimetry, soil drying, liquid infusion into porous membranes, and oxygen accumulation within polymer electrolyte membranes cueto2008nonlocal ; bazyar2018liquid ; sahimi2011flow ; berg2012stability ; bethke1991long ; dandekar2013petroleum ; benson2008carbon ; macminn2010co ; saadatpoor2010new ; bandara2011pore ; neufeld2009modelling ; bear2010modeling ; dawson1997influence ; levy2003modelling ; kueper1991two ; lee2016porous ; carmo2013comprehensive . The ability to accurately predict the displacement pathway of the non-wetting fluid is critically important in all of these cases sahimi2011flow ; rabbani2018suppressing . For instance, whether the fluid flows through a compact, stabilized front or a narrow and fingered channel impacts how much contaminant penetrates into an aquifer, how much oil can be recovered from a reservoir, and how much water remains in a dried soil cottin2011influence ; singh2003dynamic ; datta2014mobilization ; datta2013drainage .

Different displacement behaviors can arise during drainage. For a homogeneous disordered medium of uniform wettability and uncorrelated pore sizes, these behaviors are predictable using two non-dimensional parameters: the viscosity ratio $\text{M}\equiv\mu_{\text{nw}}/\mu_{\text{w}}$ and the capillary number $\text{Ca}\equiv\mu_{\text{nw}}(Q/A)/\gamma$, where $\mu_{\text{nw}}$ and $\mu_{\text{w}}$ are the non-wetting and wetting fluid viscosities, $Q$ is the fluid flow rate through a cross-sectional area $A$ of the medium, and $\gamma$ is the interfacial tension between the two fluids lenormand1988numerical ; yortsos1997phase ; xu1998invasion . Many processes are characterized by $\text{M}>1$ and $\text{Ca}\ll 1$; under these conditions, capillary forces dominate, and the non-wetting fluid cannot enter a pore of diameter $a$ until the fluid pressure reaches a threshold $\sim\gamma/a$. Therefore, the non-wetting fluid displacement proceeds one pore invasion at a time. At each time, the fluid invades the largest pore accessible to it, which is characterized by the lowest capillary pressure threshold. The flow behavior is thus determined by pore-to-pore variations in the pore size, resulting in a displacement process known as capillary fingering (CF) that is characterized by a ramified and disordered pathway lenormand1989capillary ; lenormand1983mechanisms ; lenormand1985invasion ; mayer1965mercury ; xu2008dynamics ; krummel2013visualizing ; toledo1994pore ; mason1986meniscus ; joekar2012analysis ; maaloy1992dynamics ; martys1991critical ; xu1998invasion .

Many naturally-occurring and synthetic porous media are not homogeneous, however. For example, shales, sandstones, and soils are typically heterogeneous, with smooth gradients or sharp discontinuities in pore size both along and orthogonal to the fluid flow direction ringrose1993immiscible ; schaetzl2015soils ; ashraf2019capillary . Sharp pore size stratification has been shown to alter the fluid pathway during drainage datta2013drainage ; yokoyama1981effects ; lake1981taylor ; chatzis1995investigation , yet the influence of a smooth gradient in pore sizes is still unclear. Theoretical calculations, numerical simulations, and indirect experimental evidence suggest that an applied pressure gradient can modify the fluid pathway meakin1991 ; meakin1992 ; chaouche1994 ; xu1998 , and a pore size gradient has been conjectured to play a similar role yortsos2001 , but this conjecture has not been directly tested in experiments. Recent investigations of viscous fingering, a distinct fluid displacement behavior that arises for $\text{M}<1$, demonstrate that gradients and spatial correlations in pore size can indeed strongly impact the geometry of the fluid pathway rabbani2018suppressing ; al2012control ; jackson2017stability ; pihler2012suppression ; biswas2018drying . However, how a pore size gradient impacts capillary fingering remains unknown. As a result, accurate prediction of fluid displacement pathways is still elusive for many real-world applications.

Here, we use microfluidic porous media and computational pore-network models to investigate how capillary fingering is influenced by a pore size gradient. We find that the non-wetting fluid displacement behavior depends sensitively on the direction and the magnitude of the gradient, and for a sufficiently large gradient, capillary fingering is completely suppressed. Instead, if the non-wetting fluid flows down the gradient, it propagates via a uniform, stabilized front, while if the fluid flows up the gradient, it propagates through a single, unstable, fingered channel. This behavior also depends on the relative amount of disorder in the geometry of the medium; we demonstrate that the fluid displacement can be described by a single non-dimensional parameter that quantifies the competition between a pore size gradient and pore-scale disorder. Moreover, by analyzing capillary forces at the pore scale, we develop a geometric criterion that predicts when capillary fingering is completely suppressed. Our results thereby help elucidate how diverse flow pathways can arise due to pore size gradients in disordered media.

II Experimental materials and methods

To systematically study the influence of pore size gradients and disorder on drainage, we design polydimethylsiloxane (PDMS) microfluidic devices with controlled pore geometries. The devices are comprised of a rectangular channel of width $W=$28\text{,}\mathrm{mm}\text{/} and length $L=$25\text{\,}\mathrm{mm}\text{/} containing a two-dimensional (2D) array of cylindrical pillars that span the channel height H=166~{}\pm~{}$$8\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{/}, where the uncertainty reflects experimental variations in device fabrication. Importantly, because $H\ll W$, the fluid flow is effectively 2D. We arrange the pillars on a square diagonal lattice with a diagonal length $D$ of $1.6\text{\,}\mathrm{mm}\text{/}$, as shown in Figure 1. The lattice has 31 rows with either 17 or 18 pillars per row. The pillar diameters determine the pore diameters: the diagonal distances between neighboring pillars represent the pore throats. We denote the diameters of the pores throats in row $n$ by $a_{n,i}$, where $i$ indexes the pores along the row, and the mean diameter by $\bar{a}_{n}$. Our experiments explore the range 111\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{/}$\leq\bar{a}_{n}\leq$621\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{/}, with an uncertainty of $\pm$ $1\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{/}$ due to experimental variability in device fabrication. We also introduce a triangular area at the inlet of each device to ensure that the flow is uniform before it reaches the first set of pillars rabbani2018suppressing .

We fabricate the devices using standard soft lithography. First, we design silicon molds for the devices using photolithography with a positive photoresist and deep reactive ion etching. We then cast the devices using PDMS (Dow Corning Sylgard 184), with a cross-linker to elastomer ratio of 1.5 to 10. We heat the castings for 20 minutes at $150\text{\,}\mathrm{\SIUnitSymbolCelsius}\text{/}$ and treat them for 30 seconds using a Corona treatment. Finally, we attach a glass slide to each casting and leave them to bond further overnight at $65\text{\,}\mathrm{\SIUnitSymbolCelsius}\text{/}$. The finished PDMS surfaces are hydrophobic and oleophilic, as previously described rabbani2018suppressing ; our images (e.g. in Figs. 2-3) indicate that interactions with the PDMS primarily dictate the mean curvature of the fluid-fluid interface in the experiments.

Prior to each experiment, we saturate the device with the wetting fluid, silicone oil of viscosity $\mu_{\text{w}}$ = $5\text{\,}\mathrm{mPa}\text{\,}\mathrm{s}\text{/}$. We then inject the non-wetting fluid, 76.5 vol$\%$ glycerol in water colored with food dye; the interfacial tension between the wetting and non-wetting fluids is \gamma\approx~{}$$30\text{\,}\mathrm{mN}\text{/}\mathrm{m}. The non-wetting fluid has viscosity $\mu_{\text{nw}}$ = $50\text{\,}\mathrm{mPa}\text{\,}\mathrm{s}\text{/}$; therefore, our experiments probe $\text{M}=10$. We use a Harvard Apparatus Pump 11 Elite syringe pump to impose a constant injection flow rate $Q=0.003$ mL/h, which corresponds to $\text{Ca}\equiv\mu_{\text{nw}}(Q/WH)/\gamma\approx 3\times 10^{-7}$, well in the capillary fingering regime for homogeneous disordered media. Our experiments confirm this expectation: we observe that the non-wetting fluid invades the pore space sequentially, one pore invasion at a time.

As the non-wetting fluid flows, we image the evolution of the flow pathway using a mounted digital camera taking 4K resolution images between 1 and 2 frames per minute. To characterize the structure of the resulting displacement pathway, we use the binarized images to determine the fraction of the pore space volume that is occupied by the non-wetting fluid, $V_{\text{F}}$, focusing on a field of view four rows away from the lateral boundaries and the entrance of the medium to minimize boundary effects. Moreover, to minimize exit effects, we measure $V_{\text{F}}$ when the non-wetting fluid first reaches the halfway point of the medium.

III Influence of a gradient on drainage behavior

We first test a disordered, gradient-free porous medium. To incorporate a controlled amount of disorder in the medium, we fabricate porous media with pillars whose diameters are chosen from a normal distribution. We choose the pillar diameters such that the pore diameters are also given by a normal distribution with a mean of $111\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{/}$ and with a standard deviation $\approx\sigma$ = $3.25\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{/}$, such that the difference between the maximum and minimum pore diameters equals $2\sigma$. Hence, the quantity $\sigma$ provides a measure of the disorder in the medium. Because $\text{Ca}\ll 1$, we expect the drainage to proceed via capillary fingering. Consistent with this expectation, the fluid displacement proceeds one pore invasion at a time, resulting in a ramified and disordered pathway with $V_{\text{F}}\approx$ 0.45 as shown in the rightmost panel of Fig. 1. This pathway is morphologically similar to previous observations of capillary fingering pathways, which can have $V_{\text{F}}$ ranging from $\approx 0.3$ to 0.7 lenormand1988numerical ; chen2017 .

We next test the influence of a pore size gradient on the fluid displacement. To define a gradient, we fabricate porous media with pillars whose diameters, averaged across each row, increase by a fixed amount $|g|$ per row along the imposed flow direction. As a result, the mean pore diameter decreases by $|g|$ from one row to the next, as schematized in the top row of Fig. 3 chung2017enhancing ; chaouche1994invasion ; rabbani2018suppressing . The mean pore diameter in row $n$ is thus given by $\bar{a}_{n}=\bar{a}_{1}+g(n-1)$, where the parameter $g$ quantifies the gradient; in general, a negative or positive value of $g$ indicates that the pore diameters are decreasing or increasing along the flow direction, respectively. We maintain a controlled amount of disorder in the medium by choosing the size of each pore in a given row $n$ from a normal distribution of sizes centered around $\bar{a}_{n}$ and with the difference between the maximum and minimum pore diameters again equal to $2\sigma$. To isolate the influence of the gradient, we fix $\sigma$ = $3.25\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{/}$ and test four different values of $g$: \pm$$10\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{/} and \pm$$16.7\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{/}.

Remarkably, although $\text{Ca}\ll 1$, well within the capillary fingering regime for homogeneous media, the addition of the pore size gradient completely suppresses capillary fingering. Instead, we observe two distinct flow pathways. In both cases, the non-wetting fluid displaces the wetting fluid sequentially, one pore invasion at a time; however, the macroscopic flow pathway is starkly different from the ramified and disordered pathway characteristic of capillary fingering. For $g<0$, the non-wetting fluid propagates via a uniform, stabilized front, as shown in the left panel of Fig. 2, ultimately yielding a non-wetting fluid volume fraction $V_{\text{F}}\approx 0.9$. We term this displacement behavior stable displacement. By contrast, for $g>0$, the non-wetting fluid propagates through a single, unstable, fingered channel, as shown in the right panel of Fig. 2, ultimately yielding a non-wetting fluid volume fraction $V_{\text{F}}\approx 0.1$. We term this displacement behavior unstable fingering. Interestingly, these displacement behaviors arise for both values of $|g|$, and the final $V_{\text{F}}$ appears to only depend on the sign of $g$, not its magnitude, in this range of $|g|$.

We can understand the origin of these displacement behaviors by analyzing the pore-scale capillary forces for a fluid interface trapped within a pore in a given row. Because $\text{Ca}\ll 1$, the non-wetting fluid pressure must exceed the capillary pressure threshold $\sim\gamma/a$ to enter a pore of diameter $a$. For the media with a pore size gradient, the smallest pore in row $n$, characterized by the largest capillary pressure threshold, has diameter $\bar{a}_{n}-\sigma$; by contrast, the largest pore in row $n+1$, with the smallest capillary pressure threshold, has diameter $\bar{a}_{n+1}+\sigma$. Because $\bar{a}_{n+1}-\bar{a}_{n}=-|g|$, the smallest pore in row $n$ is larger than the largest pore in row $n+1$ when $|g|/2\sigma>1$, as in the experiments, which have $|g|/2\sigma\approx 1.5$ and 2.6. Hence, when this macroscopic criterion holds—independent of the exact value of $|g|$—the pores along each row remain separated in diameter from those in the adjacent rows. Consequently, the capillary pressure threshold is lower for pores in row $n$ compared to pores in row $n+1$. The pore size gradient therefore dominates over disorder in determining the non-wetting fluid pathway: the non-wetting fluid fills all pores in a given row before proceeding to the next, ultimately leading to stable displacement as observed experimentally (Fig. 2, left).

A converse argument holds for the case of $g>0$: when $|g|/2\sigma>1$, the largest pore in a given row $n$ is smaller than the smallest pore in row $n+1$, independent of the exact value of $|g|$. The magnitudes of the corresponding capillary pressure thresholds again remain separated between adjacent rows; the pressure threshold is then higher for pores in row $n$ compared to pores in row $n+1$. The pore size gradient thus dominates over disorder in determining the non-wetting fluid pathway again: in this case, instead of laterally filling a given row, the non-wetting fluid successively fills neighboring pores in adjacent rows along the flow direction. As a result, the fluid propagates through the medium in a single thin channel approximately one pore wide, ultimately leading to unstable fingering as observed experimentally (Fig. 2, right).

IV Competition between a gradient and disorder

This geometric argument predicts that the fluid displacement pathway only depends on the sign of $g$, not its magnitude, when the non-dimensional parameter $|g|/2\sigma>1$—resulting in stable displacement or unstable fingering, as observed in the experiments. The threshold $|g|/2\sigma=1$ accordingly represents a macroscopic geometric criterion describing when capillary fingering is completely suppressed, even under conditions in which capillary fingering would typically be predicted for homogeneous disordered media. By contrast, the fluid displacement must proceed by capillary fingering in the limit that $|g|/2\sigma$ approaches zero, in which pore size disorder solely determines the displacement pathway. We thus expect a transition between these displacement behaviors as $|g|/2\sigma$ increases from zero: the pore size distributions and the corresponding capillary pressure thresholds characterizing adjacent rows overlap less and less (top row of Fig. 3, left to right). The gradient then plays an increasing role in determining the displacement pathway, with stable displacement or unstable fingering increasingly dominating for increasing values of $|g|/2\sigma$, depending on the sign of $g$, up to the threshold $|g|/2\sigma=1$.

We quantitatively test this prediction by performing experiments on devices with different values of $g$ and $\sigma$ spanning over three decades in the non-dimensional parameter $|g|/2\sigma$. Representative images showing the morphology of the fluid displacement pathway for $g<0$ are shown in the middle row of Fig. 3. Consistent with our hypothesis, we find that as the relative magnitude of the gradient increases, the non-wetting fluid displacement becomes increasingly uniform for $g<0$, ultimately leading to a stabilized front (Fig. 3, left to right). Conversely, for $g>0$, we find that as the relative magnitude of the gradient increases, the fluid displacement becomes increasingly unstable, ultimately leading to propagation through a single fingered channel.

The full set of our measurements of $V_{\text{F}}$ are shown by the large symbols in Fig. 4. Consistent with our hypothesis, $V_{\text{F}}$ increases or decreases as $|g|/2\sigma$ increases above zero for $g<0$ and $g>0$, respectively. Surprisingly, these transitions are highly sensitive to the magnitude of the gradient: the non-wetting fluid volume fraction begins to appreciably deviate from the capillary fingering value when $|g|/2\sigma$ is as small as $\sim 10^{-2}$, as shown in Fig. 4. Moreover, consistent with our hypothesis, the flow pathway eventually reaches the stable displacement and unstable fingering limits—which we measure to have $V_{\text{F}}\approx 0.9$ and $\approx 0.1$—for $g<0$ and $g>0$ respectively, as $|g|/2\sigma$ increases above one. We confirm these limiting values by also testing media with a pore size gradient ($|g|=$10\text{,}\mathrm{\SIUnitSymbolMicro m}\text{/}) but no disorder ($\sigma=0$). Additionally, to confirm that the results are insensitive to the choice of M in this regime of $\text{Ca}\ll 1$, we test different viscosity ratios $\text{M}=4$ and 40 using $65.6$ and $89.1$ vol$\%$ glycerol as the non-wetting fluid, respectively. In both cases, the gradient completely suppresses capillary fingering, as exemplified for the case of $g<0$ in Fig. [3](#S4.F3), and we measure $V_{\text{F}}=0.932$ and $0.893$ for $g=$-10\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{/} and $V_{\text{F}}=0.166$ and $0.161$ for $g=$10\text{,}\mathrm{\SIUnitSymbolMicro m}\text{/}$$, respectively. Finally, we find qualitatively similar results using different image analysis protocols as indicated by the squares and circles in Fig. 4, further indicating that our findings are robust. The slight discrepancy between the measured volume fraction in the stable displacement limit and the expected value of $V_{\text{F}}=1$ arises from experimental boundary effects. These are exemplified in Fig. 2 (Left) and in Fig. 3 (Middle Row), which show unfilled pores that occur at the boundary of the medium.

To describe the variation of $V_{\text{F}}$ with $|g|/2\sigma$, we extend the geometric argument presented in Section II. We again analyze the pore-scale capillary forces for a fluid interface trapped within a pore in row $n$, which represents the downstream tip of the non-wetting fluid interface. However, to more generally represent the pore structure, we describe the distribution of pore diameters $a_{n}$ in the row by the probability density function $p_{n}$. We denote the minimum and maximum pore diameters in row $n$ by $a_{n,\text{min}}$ and $a_{n,\text{max}}$, and the cumulative distribution function of $p_{n}$ by $\Phi_{n}$; $p_{n}(a_{n})$ then represents the probability that a pore in row $n$ has size $a_{n}$, while $\Phi_{n+1}(a_{n})$ represents the fraction of pores in row $n+1$ that are smaller than $a_{n}$. For a gradient-free medium, all the $p_{n}$ are the same and are thus independent of $n$, while as $|g|/2\sigma$ increases from zero, $p_{n}$ and $p_{n+1}$ begin to overlap less and less. We assume that the non-wetting fluid randomly samples pores as it invades the pore space, and that pore invasion only occurs at the downstream tip of the fluid interface. The non-wetting fluid volume fraction can then be approximated by the probability that a pore randomly selected from row $n$ is larger than a pore randomly and independently selected from row $n+1$:

[TABLE]

For comparison with the experiments, we represent $p_{n}$ by the normal distribution and thereby obtain a full analytic prediction for the non-wetting fluid volume fraction:

[TABLE]

Here, $\xi$ is an integration variable and thus, Equation 2 demonstrates the central role of the non-dimensional parameter $g/2\sigma$ in describing the flow dynamics. We numerically solve Eq. 2 by replacing $a_{n,\text{min}}$ and $a_{n,\text{max}}$ by $-\infty$ and $+\infty$, respectively. This solution is shown by the solid lines in Fig. 4.

The comparison between the solution to Eq. 2 and the experimental data indicates that our geometric argument qualitatively captures the transition from capillary fingering at low $|g|/2\sigma$ to stable displacement ($g<0$) or unstable fingering ($g>0$) as $|g|/2\sigma$ increases to above one. However, this argument underpredicts the sensitivity to the magnitude of the gradient: the predicted $V_{\text{F}}$ begins to appreciably deviate from the capillary fingering value only when $|g|/2\sigma\sim 10^{-1}$, an order of magnitude larger than in the experiments. This discrepancy reflects the simplifying assumptions made in this model—specifically, that the non-wetting fluid randomly samples pores as it invades the pore space, and that pore invasion only occurs at the downstream tip of the fluid interface. These assumptions neglect the sizes of and the connectivity between individual adjacent pores; they also neglect the possibility of pore invasion upstream from the downstream tip of the fluid interface, as well as wetting fluid trapping throughout the pore space, which decreases the accessibility of pores. Fully describing the influence of a pore size gradient therefore requires a more accurate description of the pore space structure and the pore-scale displacement dynamics.

V Pore-scale computational model

We hypothesize that considering pore-scale capillary forces can describe the transition to capillary fingering with increasing disorder, but only when the sizes and connectivity between individual pairs of pores in the full medium are taken into account. To test this idea, we use a computational pore-network model of invasion percolation with trapping that explicitly considers the pore-to-pore variation in capillary pressure thresholds for a specified medium geometry lenormand1985invasion ; birovljev1991gravity ; masson2016fast ; wilkinson1986percolation ; yortsos2001delineation ; xu1998invasion .

We first generate a grid of pores arranged on a square diagonal lattice with the same structure as the experimental devices. Using this grid, we use a graph object for the pore network. Pore bodies become nodes that are either “filled” by the non-wetting fluid or are “not filled”—that is, filled with the wetting fluid instead—and pore throats become edges connecting two bodies/nodes. The pore body and throat sizes are generated statistically, using specified values of $\sigma$ and $g$, in a manner similar to the experiments. Importantly, the pore throat sizes define the capillary pressure thresholds such that larger pore throats are more easily invaded than smaller ones due to the resultant differences in their capillary pressure thresholds.

At the beginning of each simulation, all the pores in the first row are filled with the non-wetting fluid. At each time step, we determine the connected component clusters of “not filled” pore bodies; the boundaries with these clusters delineate the invading non-wetting fluid interface or trapped wetting fluid regions. We then find the largest pore throat along the invading non-wetting fluid interface and fill the corresponding pore body, keeping trapped wetting fluid regions unchanged to model an incompressible fluid. The model then iterates through time steps until the non-wetting fluid reaches the outlet. Similar to the experiments, we characterize the resultant flow pathway by calculating the fraction of the pore space that is occupied by the non-wetting fluid, $V_{\text{F}}$, including both pore bodies and pore throats as in the experiments.

The non-wetting fluid pathways obtained using the pore-network simulations show excellent agreement with the experimental observations. Similar to the experimental observations described in Fig. 3, we find that for $g<0$ and increasing $|g|/2\sigma$, the non-wetting fluid displacement becomes increasingly uniform for $g<0$, ultimately leading to a stabilized front as shown in the bottom row of Fig. 3. Moreover, the flow pathway is insensitive to the choice of $\sigma$ for $|g|/2\sigma>1$. Conversely, for $g>0$, we find that as the relative magnitude of the gradient increases, the fluid displacement becomes increasingly unstable, ultimately leading to propagation through a single fingered channel.

As a final test of our hypothesis, we run pore-network simulations for geometries with different values of $g$ and $\sigma$ spanning over three decades in $|g|/2\sigma$, similar to the experiments. The mean value of $V_{\text{F}}$ obtained for each $|g|/2\sigma$ is shown by the small symbols in Fig. 4, while the full range of $V_{\text{F}}$ determined in the simulations is shown by the light shaded regions. As observed in the experiments, as $|g|/2\sigma$ increases, $V_{\text{F}}$ quickly varies from the two-dimensional capillary fingering value $\approx 0.5$, appreciably increasing ($g<0$) or decreasing ($g>0$) when $|g|/2\sigma\sim 10^{-2}$. Moreover, as observed in the experiments, the flow pathway eventually reaches the stable displacement and unstable fingering limits $V_{\text{F}}\approx 1.0$ and $\approx 0.05$ for $g<0$ and $g>0$, respectively, as $|g|/2\sigma$ increases above 1. The pore-network simulations thus establish that the macroscopic geometric criterion $|g|/2\sigma=1$ (dashed line in Fig. 4), which does not explicitly consider the pore-scale connectivity of the medium or the full pore invasion dynamics, can still describe when capillary fingering is completely suppressed.

Intriguingly, the values of $V_{\text{F}}$ obtained in the simulations are highly variable, similar to the experiments, which show a considerable amount of scatter. We do not observe any apparent correlation between $V_{\text{F}}$ and the value of $g$ or $\sigma$ used in a simulation (Appendix, Fig. A.1), confirming the governing role of the non-dimensional parameter $|g|/2\sigma$ in describing the displacement behavior. Instead, the variability in $V_{\text{F}}$ reflects the different possible pore structures that can be generated for a given value of $|g|/2\sigma$; these pore structures individually determine the different non-wetting fluid displacement pathways. Given this variability, we find good agreement between the experimental measurements and the pore-network simulations. Moreover, the variability in $V_{\text{F}}$ decreases as $|g|/2\sigma$ increases for both the experiments and the simulations: the pore size gradient increasingly dominates over disorder in determining the non-wetting fluid pathway. Together, these results support our hypothesis that for $|g|/2\sigma>1$, pore-scale capillary forces can describe the non-wetting fluid pathway, but only when the full pore space structure and the pore-scale displacement dynamics are taken into account.

VI Conclusion

Using experiments with microfluidic porous media and simulations with pore-network models, we demonstrate how capillary fingering is influenced by the competition between a pore size gradient and pore-scale disorder. We find that as the non-dimensional parameter $|g|/2\sigma$ increases, the non-wetting fluid displacement behavior quickly transitions away from the gradient-free limit of capillary fingering, even under flow conditions in which capillary fingering would typically be predicted for homogeneous disordered media. Capillary fingering is eventually completely suppressed when $|g|/2\sigma=1$, above which the non-wetting fluid flows via either stable displacement or unstable fingering, depending on the direction of the gradient. By analyzing capillary forces at the pore scale, we develop a geometric criterion to predict these limiting behaviors. Indeed, because the displacement behaviors studied here are controlled by capillary forces at the pore scale, we expect that our results are also applicable to porous media with wettability gradients. This work thus expands our understanding of drainage in porous media, which is typically described only by the two non-dimensional parameters, M and Ca lenormand1988numerical . It also suggests a way to control flow behavior in homogeneous porous media through the introduction of a pore size or wettability gradient—for example, via the controlled deposition of solute through the pore space bradford2008 .

Acknowledgements.

Correspondence should be addressed to S.S.D. at ssdatta@princeton.edu. It is a pleasure to acknowledge N. Bizmark and H. J. Cho for helpful feedback on the manuscript, and R. K. Prud’homme and H. A. Stone for stimulating discussions. We acknowledge use of the Princeton Institute for the Science and Technology of Materials (PRISM) cleanroom, and the Princeton Institute for Computational Science and Engineering for computer cluster access. This work was supported by start-up funds through Princeton University. N.B.L. was supported in part by the Mary and Randall Hack Graduate Award of the Princeton Environmental Institute. This material is also based upon work supported by the National Science Foundation Graduate Research Fellowship Program (to C.A.B.) under Grant No. DGE-1656466. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of the National Science Foundation.

Appendix A Influence of gradient and disorder on $V_{\text{F}}$

To verify that $g/2\sigma$ is an appropriate nondimensional parameter, we show that the volume fraction $V_{\text{F}}$ obtained in the simulations does not independently depend on $g$ or $\sigma$ for a given value of $g/2\sigma$. This finding is exemplified in Figure A.1: $V_{\text{F}}$ is not correlated with either $g$ (left column) or $\sigma$ (right column) for three different values of $g/2\sigma=1$, $0.1$, and $0.01$ $\pm~{}2.5\%$ (rows a, b, and c respectively). Furthermore, we calculate the Pearson correlation coefficient $\rho_{A,B}$, which is zero when variables $A$ and $B$ are uncorrelated and $\pm 1$ when they are perfectly correlated or anti-correlated. For columns a, b, and c, we find $\rho_{V_{\text{F}},g}=0.044$, $-0.009$, and $-0.003$, and $\rho_{V_{\text{F}},\sigma}=0.049$, $-0.007$, and $-0.001$, indicating a lack of correlation.