Competing Gravitational and Viscous Effects in 3D Two-Phase Flow Investigated With a Table-Top Optical Scanner

1 Introduction

Fluid dynamics in complex structures is a rapidly growing research field, considering scales ranging from the micro-scale, such as fluid transport in cells, to the macro-scale, including ocean currents and river networks. In some cases such systems enable visual study, such as monitoring of river systems through satellite imagery, but often the dynamics are hidden from view, underground or inside opaque bodies. Experimental examination of such processes is and has been an active field of study for well over a century [1]. For simplicity, both experimentally and theoretically, most of the work has been conducted in 2D systems, usually through the adaptation of some form of the Hele-Shaw cell [2]. To study 3D geometries, one can for example use X-rays [3] or nuclear magnetic resonance [4, 5], but each technique has its limitations, both in terms of what one can image and in how one can interact with the experiment. Variations on other 3D set-ups have existed at least since the 1980s [6, 7], but access to modern computers, optics, and cameras opens up for methods that previously were practically unfeasible.

In this paper we report results of experiments, in which a viscous, dense liquid is injected from above into a cubic box, filled with a synthetic porous medium, displacing a less dense, less viscous liquid. The competition between viscosity and gravity causes a transition from viscosity-stabilized flow to gravitationally unstable fingering at a critical distance from the injection point. This distance increases with the injected flow rate. Observation of the complex 3D flow patterns is enabled by a newly developed table-top experimental set-up based on optical index matching and laser-induced fluorescence.

The paper is organized as follows. In Section 2 we describe the design and configuration of the optical 3D-scanner. In Section 3 we present the first flow experiment conducted with the scanner, including methodology and theory. In Section 4 we present and analyze all results. A discussion and conclusion follow in Section 5.

2 Design and Configuration

2.1 Optical 3D Scanner Based on Index Matching

Our approach employs an underlying principle already in use in several labs [8–19], based on optical index matching and laser-excited fluorescence, known as RIM (refractive index matching methods) [20]. The advantage of our set-up is the hands-on approach, with the experimental dynamics being directly observable by eye, while remaining on a manageable scale, with the whole experiment fitting on a standard-sized table. The set-up can be assembled mostly from available products at a reasonable cost. To put this into context, our set-up costs around 10,000 EUR, which is one to two orders of magnitude less than a CT-system. However, its implementation is rather challenging, with many interconnected components that all must be controlled and accounted for. This article should therefore serve as a reference guide for anyone wishing to construct a similar set-up.

The method relies on images being captured by conventional digital cameras and the use of optical lenses, filters, and laser-induced fluorescence. Although the experiments remain within the domain of flow in porous media, the method relies heavily on optical principles. As key elements, the set-up uses two cameras to assure the image quality and increase the possible sample size, together with separate actuators for the cameras and the laser to maintain the imaging geometry and allow for constant focus at a wide aperture.

In particular, the method used here differs from methods used in other labs [8–19] in that the cameras and laser are moved by separate linear actuators, allowing for preserved imaging geometry and focus. By also using two cameras we achieve twice the possible scanning depth, allowing for larger systems. Some preliminary experimental results and a theoretical analysis are given in [21].

The refractive index of a fluid or solid, n, is the ratio of the velocity of light in vacuum to that in the medium. When media are index-matched, they are rendered transparent to one another, and the interfaces between the phases become invisible. This means that light can pass through the medium, unhindered and without distortions.

Figure 1A shows a schematic of the experimental set-up. At the center is a cubic flow cell, containing a synthetic porous medium. It is sealed by glass walls on the four sides, and by plexiglass plates with embedded valves at the center at the top and bottom. At the onset of each experiment, the sample is fully saturated by a fluid with a refractive index that as closely as possible matches that of the porous medium. To the naked eye, the flow cell appears as a homogeneous, solid block of glass, rendering the porous medium invisible. A flow experiment is initiated by injecting another, immiscible fluid with matching refractive index through the valve at the top or at the bottom. To image the flow, we introduce different fluorescent dyes in the fluid phases. The fluorescence is induced by a laser, which illuminates a vertically oriented 2D light sheet. The fluid phases appear as patches with colors corresponding to the dyes. The porous medium remains dark, as it is transparent to the light. The laser is driven in the direction normal to the illuminated sheet, using a stepper-motor driven linear actuator. The illuminated sheet is imaged, frame by frame, by the two digital cameras, which are also moved by linear actuators. The raw data thus constitute a stack of 2D images, as illustrated in Figure 1B. When organized in sequence, they provide a 3D rendering of the cell, mapping the solid and fluid phases, respectively.

FIGURE 1

FIGURE 1. (A) Experimental set-up. (B) Illustration of the unprocessed, raw data. (C) Example of a fully segmented 3D scan from a two-phase flow experiment in an artificial porous medium made from spherical glass beads.

2.2 Index-Matched Phases and Fluorescent Dyes

The experimental approach depends on the index matching of the solid and the fluid phases. However, in our configuration the refractive indices do not match perfectly, and the cell will contain some contamination. This leads to loss of light intensity and distortions that in turn reduce the distance from which the illuminated sheet can be projected. This sets limits on the image quality and the possible system size. In this particular set-up, these effects are countered by using notch filters, blocking scattered laser light, and two cameras, each of which images one half of the sample. Here it should also be noted that the pore scale comes into play, as smaller pore size will lead to a higher concentration of interfaces between the phases, which leads to more scattering, and thus to a smaller possible system size. Figure 2 demonstrates the effects of imperfect index matching. Close to the sample surface there is little distortion and loss of image quality, but at some distance the image deteriorates significantly.

FIGURE 2

FIGURE 2. Images taken at different distances into the sample. Here we see two fluid phases, red canola oil and green glycerol, in a medium of randomly packed 3 mm borosilicate glass beads. Because of the imperfect index matching, the image quality deteriorates as the illuminating light sheet is moved further into the sample.

Finding an appropriate combination of solid, immiscible fluids, and fluorescent dyes proved to be a challenge in the development of the experimental set-up. Generally, the refractive index is a function of both temperature and wavelength. The index matching can thus be tuned by controlling the temperature and by choosing a laser and dyes that emit at optimal wavelengths for the selected phases. Some fluids offer the possibility to tune the refractive index chemically, for instance by controlling the fractions in a glycerol-water mixture [22]. Tuning is more challenging with multi-phase experiments, as two or more fluid phases will not necessarily match in refractive indices with a solid phase at the same temperature. There are companies that supply tailored fluids for index-matching purposes, but they cater to microfluidics applications and their fluids are too expensive for the table-top scale of the set-up described here. The degree of index matching strongly affects the possible image quality, as well as how far into a sample it is possible to image. It therefore also influences what types of experiments are possible. In the set-up presented here, experimental simplicity is favored, leaving out temperature control and chemical tuning.

The method for deriving a good fluid pair was a mixture of literature search and trial and error (a useful reference is an online library at www.refractiveindex.info). In the end, a determining factor was the ease of handling and the safety of the fluids. For instance, silicone oil was a considered alternative, used in a similar set-up [19], but it was discarded because of the moderate toxicity and the fact that it makes it very difficult to clean the instruments between experiments.

To distinguish the porous medium and the fluid phases in the imaging, the fluids are dyed with fluorescent dyes that emit at wavelengths as widely spaced as possible, thus appearing as patches of different color. The fluids and dyes used in the experiment described here, are glycerol dyed with 10–50 mg/l Fluorescein emitting green light at 548 nm, and vegetable oil dyed with 1–10 mg/l Pyrromethene emitting red light at 650 nm, both dyes supplied by Exciton. The fluid properties are summarized in Table 1. To illuminate the cell, we use a 532 nm ZM18GF029 laser from Z-laser, pre-fitted with a prism that spreads the beam out in a 60° sheet of thickness 0.3 mm. The porous medium remains transparent to the exiting light, appearing as dark patches.

TABLE 1

TABLE 1. Reference values and experimental data for the two fluids, centered at approximately 20°C, with the span in viscosity representing ±2°C.

2.3 Imaging and Data Acquisition

The scanner connects to a computer with 128 GB of RAM and a very fast SSD disk. The mechanical operation of the linear actuators and imaging is handled by a single script. In this configuration the code is written in Python. The image data are captured through straight-forward optical imaging with standard optical lenses and cameras. In the configuration described above, there are two MER-230-168U3C 2.3MP RGB industrial vision cameras with USB-3.0 interfaces, from Get-Cameras, with 50 mm f/2.4 lenses from the same provider. The lenses are fitted with NF533-17 notch filters from Thorlabs, which filter out scattered light and reflections at wavelengths around that of the laser. These cameras have the advantage that they can be controlled directly from a script that also handles data acquisition and storage.

A full 3D scan accumulates a large amount of data. Each 3D image consists of 1,000 × 1,000 × 1,000 voxels in 8-bit 3-channel color (RGB). The raw data for each image requires about 7 GB of space that has to be either cumulatively buffered in RAM, or written to permanent storage between each scan. Therefore, the image resolution determines the possible frame rate of the scanner. We opted to configure the scanner so that the voxel resolution would be identical in all three dimensions, making a combination of the possible frame rates of the cameras and the pixel dimensions on the sensors determine the final integration time of the full scan. The cameras used can operate at up to 168 frames per second (fps), but buffering to the computer, generating a live video feed, and minimal processing limit the practical frame rate to 50 fps. A full image of 1,000 sequential 2D images takes 20 s to record. In addition, writing the image to disk and resetting the laser and camera positions bring the maximum scan rate to about one full image per minute. Although this rate is too slow for the dynamics at the onset of the faster experiments described here, it is adequate to capture the transition described in the paper. It should be noted that this rate can be increased by reducing the image resolution or buffering all data to RAM, and by recording to a hard-disk after the experiment. Because of the high viscosity of glycerol, see Table 1, leading to relatively slow dynamics, it was deemed that the scan rate was sufficiently high at 1 scan per minute. It should also be noted that the 3D images are very time consuming to process and analyze, so that gathering larger amounts of data would lead to extended processing time and require corresponding available storage space.

The imaging geometry is determined by illuminating the center cross section of the sample and shifting and refocusing one of the cameras until a desirable cropping of the sample is achieved. The other camera is then moved to an equal distance to the other side of the illuminated sheet. By measuring the pixel dimensions of an object of known size, the actual dimension of an imaged pixel can be determined.

The cameras follow the motion of the laser, but to compensate for the magnification of the illuminated plane, caused by the image being projected through a medium of higher refractive index, n, than the air between the camera lens and the sample, the laser moves at a constant speed equal to n times the speed of the cameras. In this manner a constant optical path is preserved, and the imaging magnification is kept constant. The velocities for the laser and cameras are set so that two successive 2D images correspond to focus on planes separated by one pixel dimension, i.e. they are given in Eq. 1.

$\begin{array}{cc}\hfill v_{\mathit{laser}}& =\mathit{\text{pixel dimension}}\times \mathit{\text{camera fps}}\hfill \\ \hfill v_{\mathit{cameras}}& =v_{\mathit{laser}}/n,\hfill \end{array}(1)$

where pixel dimension signifies the physical size of the individual pixels on the camera sensor. This ratio of camera speed to laser speed is an essential element of the set-up, as the constant optical path keeps the scaling of the imaging and allows the use of a wide aperture, which in turn leads to a shallow imaged depth of field, focused solely on the plane illuminated by the light sheet.

At the onset of each scan, the laser, on its carriage on the linear actuator, is outside of the sample, close to one of the vertical glass faces. The two cameras are in their corresponding positions. A scan is conducted by driving the laser sheet through the sample. The camera nearest the side the laser is in starts recording to the computer RAM, filling out a pre-initiated array. When the light sheet reaches the center of the sample, the second camera takes over. As the sheet reaches the end of the cell, the data array is written to the HD while the carriages are sent back to their original positions.

2.4 Image Processing and Segmentation

The scanner images convey an intuitive mapping of the phases to our eyes, but translating this into an automated segmentation algorithm remains challenging. As can be seen in Figure 2 such an algorithm has to deal with gradients and fluctuations in light intensity and image quality. There is no unique way to do this, but rather a whole range of possible avenues, ranging from standard image processing, to using machine learning or neural networks.

For this paper, we have used a combination of image processing with the OpenCV Python library and segmentation, analysis, and visualization with the Amira Avizo software from Thermo Fisher. The raw data consists of a 3D matrix of voxels that each have separate 8-bit RGB values. There are gradients in light intensity and image quality, as well as noise and artifacts, that have to be considered and adjusted for. The steps taken are as follows.

• Image alignment: The two cameras record the sample from opposite sides, so the recorded data from the second half needs to be flipped to match that from the first half. Although the cameras are carefully aligned, there is a possibility of slight misalignment. A transformation function, based on the findTransformECC function from OpenCV, is made by comparing the center-most images from the two cameras. This function is then applied, through the warpAffine function from OpenCV, to all the images recorded by the second camera, to rotate and stretch the image, so that the geometry is preserved throughout the sample.

• Segmentation: The data are loaded into Avizo. The segmentation contains two main steps. First, we segment out the solid phase, using a scan taken prior to the flow experiment, with the cell fully saturated with the defending phase. As the fluids contain red and green dyes respectively, this analysis is carried out solely in the image channel for the color corresponding to the defending fluid. We localize the center of mass of each glass bead and generate solid spheres with the appropriate diameter at the locations. This image of the solid matrix is then used in the segmentation of the ensuing scans. Second, the invading fluid is segmented out from the color channel corresponding to the invading fluid dye. The defending fluid is then the remaining volume. The result is a full digital mapping of the solid and the two fluid phases. Figure 3 shows an example of how the same image looks before and after segmentation.

FIGURE 3

FIGURE 3. (A) Raw data. (B) Segmented image.

3 Experiment—Competing Viscous, Capillary, and Gravitational Forces

The first experiment conducted with the set-up was a study of the interplay of viscous, capillary and gravitational forces, see Figure 4. For a quick reference, please see the video provided in the Supplementary Material section. A cubic flow cell, measuring 8 cm × 8 cm × 8 cm, was initiated with a fixed medium made up of quasi-uniform borosilicate glass spheres of 3 mm diameter supplied by Sigma Aldrich, with a reference refractive index of n ≈ 1.47, according to the supplier’s data sheet. The cell was sealed on all sides, with an inlet valve, protruding 15 mm into the matrix on top, and an outlet valve at the bottom. At the onset, the cell was fully saturated by the oil, with a red dye. We injected pure glycerol, with a green dye, at fixed flow rates from the top. The outlet was connected to a large basin, providing a quasi-fixed pressure at the outlet.

FIGURE 4

FIGURE 4. Simplified diagram of the experiment. The cell is fully saturated by liquid 2 (red) at the onset and invaded by liquid 1 (green) from the top, with a constant flow rate, Q. The outlet leads to a large-surface reservoir open to the atmosphere, in order to maintain a near constant pressure at the outlet. According to our approximate theory, the invasion starts as a sphere, R₁, until the front hits the roof, where the front gradually changes to a downward-pointing hemisphere, R₂.

We conducted experiments at different flow rates, and analyzed the obtained invasion patterns. The properties of the fluid phases are summarized in Table 1.

3.1 Theoretical Context

As a direct parallel to experiments considered in 2D systems [23–26], our theoretical starting point is that the ensuing dynamics are mainly governed by an interplay between forces dictated by the fluid properties, the wetting conditions, the geometric configuration of the matrix, gravity, and the invasion method. Depending on these parameters, the invading fluid displaces the defending fluid on a scale ranging from full displacement leading to near total saturation, to unstable branched displacement leading to low saturation. The invasion body assumes a geometrical shape which retains characteristics that can be quantified by saturation and topological features. This geometrical structure is the main result delivered by the scanner. We distinguish between the high saturation scenario, with an invasion front preserving the symmetries imposed by the injection geometry and gravity, which we call stable displacement, or low saturation, spontaneous symmetry breaking of the front and fingering, unstable displacement. For instance, if there is a density contrast between the fluids, the more dense fluid front will be gravity stabilized when invading upwards, and destabilized when invading downward. A wider spread in pore-throat size distribution will lead to larger fluctuations in the front, as the flow will favor the wide passages and leave behind pores accessible only through narrow pore throats.

In two-phase flow in a porous medium, the condition that governs if a given pore is invaded or not is whether the pressure difference between the two sides of the interface exceeds the capillary pressure threshold set up by the surface tension, γ, between the fluid phases in the local pore geometry.

For single-phase flow in a porous medium Darcy’s law on a differential form, may be written as

$\vec{v}=-\frac{\kappa }{\mu }\left(\nabla p-\rho \vec{g}\right),(2)$

where $\vec{v}$ is the filtration velocity, ρ is the fluid density, $\vec{g}$ is the gravitational acceleration, μ is the fluid viscosity, ∇p is the pressure gradient, and κ is the permeability of the medium.

For a single fluid phase entering with a constant flow rate, Q, at depth z = 0 into a vertical pipe containing a porous medium with cross section A, Darcy’s can be rewritten on a integral form as

$p\left(z\right)=p\left(0\right)+\left(\rho g-\frac{\mu Q}{\kappa A}\right)z,(3)$

where p(z) is the pressure at depth z. Here $\vec{z}$ is pointing downward, parallel with $\vec{g}$. This gives a local pressure, at depth z, from gravity, as well as a viscous contribution.

In this experiment, invasion occurs from above with a more viscous, more dense fluid. In this situation, gravity is expected to destabilize the invasion, in competition with the stabilizing viscous pressure drop [27]. In the case of high local fluid velocity, we expect stabilization of the front [23–25]. Due to mass conservation, the fluid velocity decays as a function of distance to the inlet. Based on this, we derive a theory in which we assume a stable invasion near the narrow inlet, caused by the high local velocity. The point injection with gravitational field leads to a cylindrical symmetry of the problem around a vertical axis aligned with the injection point, far enough from the side boundaries. Close to the inlet where velocity is high, viscous pressure gradients dominate over gravitational pressure gradients, and the velocity flow is expected to be symmetric under rotation around the injection point. Indeed, from experimental observations, we find it reasonable to approximate the shape of the invading fluid during this stable phase as a downward-facing hemisphere centered on the inlet. See Figure 4.

We note that the wetting condition is a natural feature to consider in this context. We conducted measurements of fluid-fluid interface tension and fluid-solid wetting angles with a drop shape analyzer and found the wetting conditions to be indeterminably neutral. We observe that the invading glycerol appears to be the non-wetting phase when it is the invading phase, going through the pore-throats with a convex front, whereas the oil also assumes a convex front shape, when it is used as the invading phase. The wetting conditions have therefore not been included in our theory.

3.2 Permeability Measurements

The permeability in the experimental context is a complex issue. Neither the traditional one-phase permeability, as defined by Darcy in Eq. 3 nor relative multi-phase permeability as it is traditionally defined [28] applies fully to this case, as the fluid saturations evolve and the fluid interface configurations are not in a steady state. To get a best approximation we have made a set-up, quite along the lines of Darcy’s experiment, as illustrated and explained in Figure 5.

FIGURE 5

FIGURE 5. Set-up for permeability measurements. The cell is open at the top, with an overflow to maintain near-constant pressure. Fluid is injected through an inlet at the bottom, with constant flow rate maintained by a syringe pump. Holes with nipples are made over two fixed intervals for preassure measurements, one with the nipple at the cell wall, and one with the nipple inserted a couple of centimeters into the matrix, to verify that the pressure gradient is the same inside the matrix, as along the walls. Pressure differences are read off over the intervals as hydrostatic fluid columns in transparent tubes, attached vertically over a ruler to measure height differences, Δh, as seen on the right. The viscous pressure gradients then can be calculated as |∇p| = ρ₁gΔh/ΔH.

Two sets of measurements were made, one for single-phase flow, and one for two-phase flow. For the single-phase flow, the cell was pre-saturated with glycerol, and the same fluid was injected at a fixed flow rate, Q at the bottom. For the two-phase flow, the cell was pre-saturated with the oil, after which glycerol was injected until the read-off tubes were all in direct contact with the invading glycerol. Any remaining oil in the tubes was removed with a syringe. Two-phase measurements were then made with glycerol injected at a fixed flow rate at the bottom. The permeability was calculated from Eq. 3 as

$\kappa =\frac{v{\mu }_1}{|\nabla p|}=\frac{Q{\mu }_1\mathrm{\Delta }H}{A{\rho }_{1}g\mathrm{\Delta }h}.(4)$

The results are presented in Table 2. As a reference, the Kozeny-Carman relation [19] puts the permeability at $\kappa \approx \frac{r^2}{45}\frac{{\varphi }^3}{{(1-\varphi )}^2}=1.1\times 10^{-8}$ m², where r = 1.5 mm is the bead radius, and ϕ ≈ 0.42 is the porosity measured from the experiments.

TABLE 2

TABLE 2. Results from single-phase and two-phase permeability measurements. The error is evaluated as the standard deviation in the measurements.

3.3 The Critical Radius

In a real porous medium the front will never be a perfect hemisphere, but will have fluctuations due to the variations in pore sizes with corresponding capillary pressure-threshold fluctuations. We therefore consider a small perturbation of the hemisphere, corresponding to an interface in a pore located a distance a in front of the hemisphere. We consider the two points, z = R on the hemisphere and a perturbation at z = R + a just in front of that point, representing a pore that has just been invaded. We compare the capillary pressure in this pore, with that at a point on the hemisphere at z = R. As an approximation using Eq. 2 in polar coordinates for the hemisphere, the change of pressure in fluid 1 and fluid 2, respectively, can be expressed as

$\mathrm{\Delta }{p}_1=p_1\left(z=R+a\right)-p_1\left(z=R\right)=\left({\rho }_{1}g-\frac{Q{\mu }_1}{A\left(R\right){\kappa }_{TP}}\right)a,(5)$

$\mathrm{\Delta }{p}_2=p_2\left(z=R+a\right)-p_2\left(z=R\right)=\left({\rho }_{2}g-\frac{Q{\mu }_2}{A\left(R\right){\kappa }_{SP}}\right)a,(6)$

where A(R) = 2πR² is the area of the hemisphere. Let p^c = p₁ − p₂ be the capillary pressure across the interface between fluid 1 and fluid 2 at a given position. The increase in capillary pressure is then

$\begin{array}{cc}\hfill \mathrm{\Delta }{p}^c& =p^c\left(z=R+a\right)-p^c\left(z=R\right)\hfill \\ \hfill & =\mathrm{\Delta }{p}_1-\mathrm{\Delta }{p}_2\hfill \\ \hfill & =\left[\left({\rho }_1-{\rho }_2\right)g-\left(\frac{Q{\mu }_1}{A\left(R\right){\kappa }_{TP}}-\frac{Q{\mu }_2}{A\left(R\right){\kappa }_{SP}}\right)\right]a.\hfill \end{array}(7)$

The capillary pressure at the perturbation will have an increase due to the hydrostatic pressure and a decrease due the viscous pressure corresponding to the first and the second term in the equation above, respectively. If Δp^c > 0, the gravitational term dominates and the invasion is more likely to continue at the perturbation, than at the lagging front. The front is thus rendered unstable. Conversely, if Δp^c < 0, the viscous term dominates and the front remains stable. We therefore use Δp^c = 0, with the front area approximated as $A\approx 2\pi R_c^2$, to predict the characteristic length scale R_c for the transition between the viscous stable and the gravitationally unstable displacement.

$R_c=\sqrt{\left[\frac{{\mu }_1}{{\kappa }_{TP}}-\frac{{\mu }_2}{{\kappa }_{SP}}\right]\frac{Q}{2\pi \mathrm{\Delta }\rho g}}.(8)$

We note that since μ₁ ≫ μ₂, the $\frac{{\mu }_2}{{\kappa }_{SP}}$ term barely contributes to the estimated values for R_c.

The theory described above predicts a length scale for the transition from stable, viscosity-dominated invasion to unstable, gravity-driven invasion. However, it does not include the role of the spatial fluctuations in the porous medium. We have described these effects in greater detail in [21]. As the invasion front progresses, the porous medium presents a distribution of capillary pressure thresholds due to the pinned menisci in the pore throats. As long as − Δp^c is negative, with a magnitude much larger than the width of the capillary pressure distribution W_t nearly every pore will be invaded and the interface is stable. The width of stable invasion fronts in quasi 2D experiments has been found to depend on the generalised fluctuation number F = −Δp^c/W_t with a transition to unstable fronts when F = 0 [23, 29–32]. In our experiments F = 0 corresponds to the situation where the invasion structure reaches the critical radius where the local fluctuations are felt, and the wider pores are preferred over narrower ones. This leads to fingering and unstable displacement, such as investigated in 2D [33].

4 Results

Figure 6 shows a 3D rendering from a scanner image. The scanner images allow for analyses from the pore scale up to the meso scale of the full flow cell, with the limits set by fine tuning of the refractive indices and a bigger system.

FIGURE 6

FIGURE 6. 3D image from the scanner, from an experiment with Q = 0.8 ml/min, captured at the moment the invading fluid reached the bottom of the cell, after 31 min. The porous medium is cut in half, with the right side removed for visual access to the invasion body.

4.1 Flow Experiments

We conducted fourteen runs, in which we increased the flow rate from 0.1 ml/min to 8.0 ml/min. Each experiment was conducted after resetting the flow cell, with a new bead configuration for the porous medium. Based on our assumptions, the stabilizing viscous pressure drops radially from the inlet in the invading fluid. This should lead to the flow being unstable, dominated by gravitational forces beyond the critical radius, R_c [5, 30, 32–34]. Figure 7 shows density projections of the 3D images from four of the experiments. As can be seen, the structures become denser and appear to assume a hemispherical shape as the flow rate is increased.

FIGURE 7

FIGURE 7. Density projections from the scanner images of the invasion body for four different flow rates. All the images are taken at the moment the fluid reaches either the bottom or the wall of the flow cell. The top row is seen from above and the bottom row is seen from the side. The projections are made by adding the layers of the 3D image, mapping each layer with ones for the invading fluid, and zeros elsewhere, into the image plane. As can be seen, the body becomes denser, with higher global saturation, as the flow rate is increased. The red crosses indicate the location of the inlet.

To identify a critical radius, R_c, as in Eq. 8, we consider the saturation near the inlet. In Figure 8 we show a plot of the saturation inside downward-facing hemispherical shells, centered at the inlet. We measure R_c as the radius where the saturation in the shell drops below 0.95. These plots support the assumption that the viscous stabilization only becomes significant beyond a given flow rate.

FIGURE 8

FIGURE 8. Plots of the saturation inside a downward-facing hemisphere with center at the inlet, imposed over the invasion structures. The points indicate the saturation of the invasion body inside 1-mm incremental hemispherical shells. The black stars are the measured transition radius, R_c, corresponding to a set threshold of saturation, S = 0.95. The black dots are values corresponding to S = 0.98 and S = 0.92, which we use as the error range for R_c.

Figure 9 shows 3D renderings of the invading structures from scans made after approximately 23 ml of invading fluid was injected, for two different flow rates, 1 and 8 ml/min. As expected, the structures become more dense and occupy a higher saturation for the higher flow rate. The higher the flow rate, the more reasonable the hemisphere assumption becomes.

FIGURE 9

FIGURE 9. 3D renderings of the scanner images of the invasion body for two different flow rates with approximately 23 ml injected in both cases. The green signifies the outline of the full body, and the red inside is segmented out as the stable part of the same body. The segmentation procedure, known as morphological opening (erosion-dilation), is a common image-processing procedure, for instance described in [35].

Table 3 shows data obtained from images of the experiments. As can be seen, the three lowest flow rates lead to very similar global saturation, while beyond a flow rate of 0.4 ml/min we see a more systematic trend toward the invasion being more stabilized, resulting in higher saturation.

TABLE 3

TABLE 3. Experimental data. For Q = 0.1 ml/min to Q = 2.0 ml/min, the corresponding images are taken the moment the invading fluid reaches the bottom of the cell. For Q = 4.0 and 8.0 ml/min, the structures reached the cell walls before reaching the bottom and were imaged at that time. Thus, data for time and injected volume are not strictly comparable. It should also be noted that the times and corresponding volumes for break-through are approximate, as scans where conducted at 2-min intervals. R_c measurements are calculated as explained in Figure 8, with the error as the larger of the difference in R_c, from the reference 0.95 saturation threshold to the 0.92 and 0.98 thresholds.

Figure 10 shows a plot of the measured transition radii vs. the flow rate, together with the theoretical prediction, Eq. 8.

FIGURE 10

FIGURE 10. The measured transition radii for the 14 flow rates, and the theoretical prediction, Eq. 8. The red dots are measurements where the radii are measured below one pore radius and the two fastest experiments, where the cell was not large enough to determine if the transition had occurred.

5 Discussion and Conclusion

The measured critical radii, R_c, fit well with the described theory, but there are not enough experiments, and the considered system is not large enough to conclude firmly and generate statistics over several decades. However, our results are significant in charting the way for experimental work on flow in 3D porus media. These experiments are inherently difficult, and there will always be limitations.

In this paper we have described the steps for the construction of a table-top set-up for imaging multiphase flows in 3D porous media. This set-up is considerably less expensive than traditional methods, such as X-ray or magnetic resonance set-ups and allows for direct optical visualization of flows without any need for image inversion algorithms. We have employed the set-up to show how viscous forces can stabilize an invasion front that would be rendered unstable due to an unfavorable density contrast (heavier fluid injected on top of a lighter one). Our set-up allowed for the complete 3D characterization of the invading structure. By considering the balance between the (stabilizing) viscous and (destabilizing) gravitational forces, we have derived an expression for flow-rate dependence of the critical radius, beyond which the front becomes unstable. Although the experimental dataset does not extend far enough for us to accurately measure the dependence of this critical radius on the flow rate, our experiments are in reasonable agreement with the theoretical prediction.

The experimental set-up illustrated here can in principle be further developed to study flows of more fluid phases, the limitation for this being the constraint of matching the index of refraction of the additional fluid phases to the ones already present. As mentioned previously, the limitations due to the imperfect matching of the refraction index can lead to a reduction of image quality, and one must take into account the desired pore-scale resolution when deciding how accurate the index matching needs to be. It should be noted that better index matching should possible. This would allow for larger system sizes, either through using a larger cell or with smaller glass beads, but we deemed the current state sufficient to start conducting experiments. We have conducted some tests with 1 mm beads, and although the possible real dimensions of the system goes down, a larger number of pores can be imaged, but as the laser sheet has a thickness of 0.3 mm and optical constraints, the pore scale detail rendering is much poorer, setting a limit on the smallest pore size that it is possible to resolve.

Keeping the optical path constant is a critical component of the set-up we describe. This is achieved by tuning the speeds of the laser sheet and the cameras accordingly, as illustrated in Eq. 1. In spite of being a relatively simple step to execute, our experience with the set-up showed early on that this key element has a dramatic effect on the overall quality of the 2D images (and consequently on the quality of the 3D reconstruction).

Our optical scanner adds to the tools of imaging capabilities available to the scientific community. In spite of its relative simplicity and low cost, it does not invalidate other techniques traditionally employed. Each imaging technique has its own applicability range and associated limitations or drawbacks. Here, the advantage of not needing radioactive materials or large facilities is counteracted by the limitation of needing fluids with matching refraction indices. In applications in which this requirement can be met, we believe this technique can be useful for a wide range of systems.

We gratefully acknowledge the support from the University of Oslo and the Research Council of Norway through its Centre of Excellence funding scheme with project number 262644.

Data Availability Statement

The raw data supporting the conclusion of this article will be made available by the authors, without undue reservation.

Author Contributions

JFB made the experimental set-up, conducted the experiments and wrote the bulk of the article. PAR partook in discussions throughout, observed some of the experiments, edited the text, and wrote contributing paragraphs. MM partook in discussions throughout, observed some of the experiments, edited the text, and wrote contributing paragraphs. RT partook in discussions throughout and edited the text. KJM participated in planning and initiating the set-up and experiments, partook in discussions throughout, observed some of the experiments, edited the text, and wrote contributing paragraphs.

Funding

The Research Council of Norway through its Centre of Excellence funding scheme with project number 262644. MM acknowledges the support from the Research Council of Norway through its Researcher Project for Young Talents funding scheme with project number 324555.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We gratefully acknowledge the suppurt of the University of Oslo and thank Mihailo Jankov for discussions and technical support.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphy.2022.936915/full#supplementary-material

References

1. Darcy H. Les fontaines publiques de la ville de Dijon: exposition et application des principes à suivre et des formules à employer dans les questions de distribution d'eau. Paris: Victor Dalmont, Libraire des Corps imperiaux des ponts et chaussees et des mines (1856).

Google Scholar

2. Hele-Shaw HS. The Flow of Water. Nature (1898) 58:34–6. doi:10.1038/058034a0

CrossRef Full Text | Google Scholar

3. Berg S, Ott H, Klapp SA, Schwing A, Neiteler R, Brussee N, et al. Real-time 3d Imaging of Haines Jumps in Porous media Flow. Proc Natl Acad Sci U.S.A (2013) 110:3755–9. doi:10.1073/pnas.1221373110

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Allen SG, Stephenson PCL, Strange JH. Morphology of Porous media Studied by Nuclear Magnetic Resonance. J Chem Phys (1997) 106:7802–9. doi:10.1063/1.473780

CrossRef Full Text | Google Scholar

5. Yan J, Luo X, Wang W, Toussaint R, Schmittbuhl J, Vasseur G, et al. An Experimental Study of Secondary Oil Migration in a Three-Dimensional Tilted Porous Medium. AAPG Bulletin (2012) 96:773–88. doi:10.1306/09091110140

CrossRef Full Text | Google Scholar

6. Clément E, Baudet C, Hulin J. Multiple Scale Structure of Non Wetting Fluid Invasion Fronts in 3d Model Porous media. J de Physique Lettres (1985) 46:1163–71.

Google Scholar

7. Frette V, Måløy KJ, Boger F, Feder J, Jøssang T, Meakin P. Diffusion-limited-aggregation-like Displacement Structures in a Three-Dimensional Porous Medium. Phys Rev A (1990) 42:3432–7. doi:10.1103/PhysRevA.42.3432

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Harshani HMD, Galindo-Torres SA, Scheuermann A, Muhlhaus HB. Experimental Study of Porous Media Flow Using Hydro-Gel Beads and LED Based PIV. Meas Sci Technol (2017) 28:015902. doi:10.1088/1361-6501/28/1/015902

CrossRef Full Text | Google Scholar

9. Moroni M, Kleinfelter N, Cushman JH. Analysis of Dispersion in Porous media via Matched-index Particle Tracking Velocimetry Experiments. Adv Water Resour (2007) 30:1–15. doi:10.1016/j.advwatres.2006.02.005

CrossRef Full Text | Google Scholar

10. Datta SS, Dupin J-B, Weitz DA. Fluid Breakup during Simultaneous Two-phase Flow through a Three-Dimensional Porous Medium. Phys Fluids (2014) 26:062004. doi:10.1063/1.4884955

CrossRef Full Text | Google Scholar

11. Stöhr M, Roth K, Jähne B. Measurement of 3D Pore-Scale Flow in Index-Matched Porous Media. Exp Fluids (2003) 35:159–66. doi:10.1007/s00348-003-0641-x

CrossRef Full Text | Google Scholar

12. Roth EJ, Mont-Eton ME, Gilbert B, Lei TC, Mays DC. Measurement of Colloidal Phenomena during Flow through Refractive Index-Matched Porous Media. Rev Scientific Instr (2015) 86:113103–9. doi:10.1063/1.4935576

CrossRef Full Text | Google Scholar

13. do Nascimento DF, Vimieiro Junior JR, Paciornik S, Carvalho MS. Pore Scale Visualization of Drainage in 3D Porous Media by Confocal Microscopy. Sci Rep (2019) 9:1–11. doi:10.1038/s41598-019-48803-z

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Ovdat H, Berkowitz B. Pore-scale Study of Drainage Displacement under Combined Capillary and Gravity Effects in Index-Matched Porous Media. Water Resour Res (2006) 42. doi:10.1029/2005wr004553

CrossRef Full Text | Google Scholar

15. Sharma P, Aswathi P, Sane A, Ghosh S, Bhattacharya S. Three-dimensional Real-Time Imaging of Bi-phasic Flow through Porous Media. Rev Scientific Instr (2011) 82:113704. doi:10.1063/1.3658822

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Kong X-Z, Holzner M, Stauffer F, Kinzelbach W. Time-resolved 3d Visualization of Air Injection in a Liquid-Saturated Refractive-Index-matched Porous Medium. Exp Fluids (2011) 50:1659–70. doi:10.1007/s00348-010-1018-6

CrossRef Full Text | Google Scholar

17. Heyman J, Turuban R, Lester DR, Méheust Y, Le Borgne T. Transverse Mixing in 3d Porous Flows: Experimental Results & Stochastic Modelling. In: Computer Methods in Water Resources (CMWR) 2018. Nantes, France (2018).

Google Scholar

18. Kang JH, Lee K-J, Nam JH, Kim C-J, Park HS, Lee S, et al. Visualization of Invasion-Percolation Drainage Process in Porous Media Using Density-Matched Immiscible Fluids and Refractive-Index-Matched Solid Structures. J Power Sourc (2010) 195:2608–12. doi:10.1016/j.jpowsour.2009.11.087

CrossRef Full Text | Google Scholar

19. Dalbe MJ, Juanes R. Morphodynamics of Fluid-Fluid Displacement in Three-Dimensional Deformable Granular Media. Phys Rev Appl (2018) 9:024028. doi:10.1103/physrevapplied.9.024028

CrossRef Full Text | Google Scholar

20. Wiederseiner S, Andreini N, Epely-Chauvin G, Ancey C. Refractive-index and Density Matching in Concentrated Particle Suspensions: A Review. Exp Fluids (2011) 50:1183–206. doi:10.1007/s00348-010-0996-8

CrossRef Full Text | Google Scholar

21. Brodin JF, Moura M, Toussaint R, Måløy KJ, Rikvold PA. Visualization by Optical Fluorescence of Two-phase Flow in a Three-Dimensional Porous Medium. J Phys Conf Ser (2022) 2241:012004. doi:10.1088/1742-6596/2241/1/012004

CrossRef Full Text | Google Scholar

22. Heyman J, Lester DR, Turuban R, Méheust Y, Le Borgne T. Stretching and Folding Sustain Microscale Chemical Gradients in Porous Media. Proc Natl Acad Sci U.S.A (2020) 117:13359–65. doi:10.1073/pnas.2002858117

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Birovljev A, Furuberg L, Feder J, Jøssang T, Måløy KJ, Aharony A. Gravity Invasion Percolation in Two Dimensions: Experiment and Simulation. Phys Rev Lett (1991) 67:584–7. doi:10.1103/physrevlett.67.584

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Frette V, Feder J, Jøssang T, Meakin P. Buoyancy-driven Fluid Migration in Porous Media. Phys Rev Lett (1992) 68:3164–7. doi:10.1103/PhysRevLett.68.3164

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Løvoll G, Méheust Y, Måløy KJ, Aker E, Schmittbuhl J. Competition of Gravity, Capillary and Viscous Forces during Drainage in a Two-Dimensional Porous Medium, a Pore Scale Study. Energy (2005) 30:861–72. doi:10.1016/j.energy.2004.03.100

CrossRef Full Text | Google Scholar

26. Toussaint R, Måløy KJ, Méheust Y, Løvoll G, Jankov M, Schäfer G, et al. Two-phase Flow: Structure, Upscaling, and Consequences for Macroscopic Transport Properties. Vadose Zone J (2012) 11:vzj2011.0123. doi:10.2136/vzj2011.0123

CrossRef Full Text | Google Scholar

27. Sabet N, Hassanzadeh H, De Wit A, Abedi J. Scalings of Rayleigh-Taylor Instability at Large Viscosity Contrasts in Porous Media. Phys Rev Lett (2021) 126:094501. doi:10.1103/PhysRevLett.126.094501

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Johnson EF, Bossler DP, Bossler VON. Calculation of Relative Permeability from Displacement Experiments. Trans AIME (1959) 216:370–2. doi:10.2118/1023-g

CrossRef Full Text | Google Scholar

29. Wilkinson D. Percolation Model of Immiscible Displacement in the Presence of Buoyancy Forces. Phys Rev A (1984) 30:520–31. doi:10.1103/physreva.30.520

CrossRef Full Text | Google Scholar

30. Méheust Y, Løvoll G, Måløy KJ, Schmittbuhl J. Interface Scaling in a Two-Dimensional Porous Medium under Combined Viscous, Gravity, and Capillary Effects. Phys Rev E (2002) 66:051603. doi:10.1103/PhysRevE.66.051603

CrossRef Full Text | Google Scholar

31. Toussaint R, Løvoll G, Méheust Y, Måløy KJ, Schmittbuhl J. Influence of Pore-Scale Disorder on Viscous Fingering during Drainage. Europhys Lett (2005) 71:583–9. doi:10.1209/epl/i2005-10136-9

CrossRef Full Text | Google Scholar

32. Måløy KJ, Moura M, Hansen A, Flekkøy EG, Toussaint R. Burst Dynamics, Upscaling and Dissipation of Slow Drainage in Porous Media. Front Phys (2021) 718. doi:10.3389/fphy.2021.796019

CrossRef Full Text | Google Scholar

33. Vasseur G, Luo X, Yan J, Loggia D, Toussaint R, Schmittbuhl J. Flow Regime Associated with Vertical Secondary Migration. Mar Pet Geology (2013) 45:150–8. doi:10.1016/j.marpetgeo.2013.04.020

CrossRef Full Text | Google Scholar

34. Auradou H, Måløy KJ, Schmittbuhl J, Hansen A, Bideau D. Competition between Correlated Buoyancy and Uncorrelated Capillary Effects during Drainage. Phys Rev E (1999) 60:7224–34. doi:10.1103/physreve.60.7224

CrossRef Full Text | Google Scholar

35. Said KAM, Jambek AB, Sulaiman N. A Study of Image Processing Using Morphological Opening and Closing Processes. Int J Control Theor Appl (2016) 9:15–21.

Google Scholar