The Imaging of Dynamic Multiphase Fluid Flow Using Synchrotron-Based X-ray Microtomography at Reservoir Conditions

Fast synchrotron-based X-ray microtomography was used to image the injection of super-critical CO2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {CO}_{2}$$\end{document} under subsurface conditions into a brine-saturated carbonate sample at the pore-scale with a voxel size of 3.64μm\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$3.64\,\upmu \hbox {m}$$\end{document} and a temporal resolution of 45 s. Capillary pressure was measured from the images by finding the curvature of terminal menisci of both connected and disconnected CO2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {CO}_{2}$$\end{document} clusters. We provide an analysis of three individual dynamic drainage events at elevated temperatures and pressures on the tens of seconds timescale, showing non-local interface recession due to capillary pressure change, and both local and distal (non-local) snap-off. The measured capillary pressure change is not sufficient to explain snap-off in this system, as the disconnected CO2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {CO}_{2}$$\end{document} has a much lower capillary pressure than the connected CO2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {CO}_{2}$$\end{document} both before and after the event. Disconnected regions instead preserve extremely low dynamic capillary pressures generated during the event. Snap-off due to these dynamic effects is not only controlled by the pore topography and throat radius, but also by the local fluid arrangement. Whereas disconnected fluid configurations produced by local snap-off were rapidly reconnected with the connected CO2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\hbox {CO}_{2}$$\end{document} region, distal snap-off produced much more long-lasting fluid configurations, showing that dynamic forces can have a persistent impact on the pattern and sequence of drainage events.

of such systems have been limited to the core-scale. However, in the last 25 years X-ray microtomography has allowed for direct imaging of the pore-space of rocks, providing input for pore-scale models and allowing for direct experimental investigation of fluid distributions Cnudde and Boone 2013;Wildenschild and Sheppard 2013).
Early results relied on either low-intensity laboratory-based or synchrotron X-ray sources, with scans typically taking several hours to acquire. Flow processes at the pore-scale tend to occur on much shorter timescales, meaning that experiments were limited to studying static flow end points, after displacement processes have finished. In recent years, however, highintensity synchrotron light sources have allowed for imaging with a high temporal resolution, sufficient for the imaging of dynamic displacement processes, with pioneering papers by Berg et al. (2013) and Armstrong et al. (2014) examining Haines jumps and capillary desaturation, respectively.
Multiphase flow is controlled by multiple different factors, including properties of the porous medium, such as the chemical composition of the system and pore topography, topology and distribution. The behaviour is also determined by the fluid properties, such as contact angle, density, interfacial tension and viscosity. The properties of the pore-space and the fluids can even be interdependent, as demonstrated by the influence that changes in surface roughness have on the contact angle (Wenzel 1936;Cassie and Baxter 1944).
As the fluid properties are dependent on system temperature and pressure (Espinoza and Santamarina 2010;Li et al. 2012;Spiteri et al. 2008;Wang and Gupta 1995), the examination of such a complex system at realistic conditions is highly attractive. Recent studies have shown that, as X-ray microtomography is non-invasive, it is possible to conduct experiments at representative subsurface conditions to investigate capillary trapping (Andrew et al. 2013(Andrew et al. , 2014cChaudhary et al. 2013), measure contact angle (Andrew et al. 2014b) and examine multiple points of changing capillary pressure during drainage (Herring et al. 2014). However, the examination of time-resolved dynamic fluid flow on a pore-by-pore level has hitherto been limited to ambient conditions Armstrong et al. 2014;Armstrong and Berg 2013;Rücker et al. 2015).
Drainage, or the injection of non-wetting phase (CO 2 ) into a porous medium saturated with wetting phase (brine), is a complex process in reservoir rocks. A qualitative understanding can be found using network modelling, where the geometry of pore-space is simplified as a lattice of pores connected by throats. Using this description, drainage can be described using invasion percolation theory where the non-wetting phase (CO 2 ) displaces the wetting phase (brine) from each pore in a sequence of pore-scale events called Haines jumps (Wilkinson and Willemsen 1983).
Haines jumps are typically assumed to be quasi-static; they are an irreversible change between two states at capillary equilibrium. The Haines jumps therefore represent the redistribution of surface energy from an old energetically unfavourable state to a new, more energetically favourable state. Each of the states therefore represents a local surface energy minimum. An event occurs when the old state ceases to be an energy minimum, and the fluid redistributes to find minimum new position of equilibrium. This allows for the sequence of events to be unambiguously predicted if the initial conditions, pore-space geometry and topology (connectedness) are known. A description of the dynamics of the transition between the two states is, however, lacking, and any change in the capillary pressure associated with the jump is assumed to be confined to the pores and throats immediately involved with the event (Blunt et al. 2002). The presence of cooperative pore-filling and non-local changes in capillarity have been demonstrated using both rate-controlled mercury injection (Yuan 1990; Mohanty et al. 1987) and the rapid 2D imaging of micromodel drainage (Armstrong and Berg 2013), although this work does not directly challenge the assumption of a series of equilibrium states.
The main goal of our study is to identify individual dynamic drainage events at the porescale, at reservoir conditions, to advance our understanding of these phenomena, and in particular to test some of the assumptions of local capillarity and quasi-static flow used in the formulation of pore-scale network (Valvatne and Blunt 2004;Øren et al. 1998) and direct simulation models (Raeini et al. 2014a). Specifically, we observe persistent fluid configurations of trapped non-wetting phase, which cannot be explained by global capillary equilibrium concepts alone.
We first outline our method by which the dynamics of multiphase fluid flow at elevated temperatures and pressures has been examined using 4D synchrotron X-ray microtomography. We focus on studying both local and non-local dynamic pore-scale events during drainage in which the injection of a non-wetting phase (supercritical CO 2 ) into brine (wetting phase) is performed in a carbonate rock. We then demonstrate how the experimental method developed in this work can be used to describe and explain in detail three pore-scale phenomena during drainage: interface recession due to capillary pressure change, Roof snap-off (Roof 1970) and distal (non-local) snap-off. The relative importance of these phenomena to the overall progress of the drainage process is beyond the scope of this study and should be the target for future research.
The terms "event" and "jump" are used synonymously to refer to the rapid changes in fluid distribution imaged in this study. The term "image" is used to refer to a reconstructed 3D volume.

Materials and Methods
The experimental apparatus is shown in Fig. 1. All imaging was performed at the Diamond Imaging Beamline I13 (Diamond Light Source, UK). The experiments were conducted on Ketton Oolite, a quarry limestone from the Upper Lincolnshire Limestone Member with a porosity of 0.2337 and a gas permeability of 2.81 × 10 −12 m 2 (2.77 D) (Weatherford Laboratories, East Grinstead, UK).
Basic information about the experiment can be found in Table 1.

Flow Boundary Conditions
To control the extremely low flow rates required for dynamic imaging, a novel flow boundary condition was used. Instead of using constant flow and constant pressure at the inlet and outlet of the core respectively, a very low-permeability porous plate was incorporated in series with the core and a constant pressure drop was applied across it. For this a hydrophilic-modified semipermeable disc (aluminium silicate, Weatherford Laboratories, Stavanger, Norway) with a permeability of 1.4 × 10 −17 m 2 (14 × 10 −6 D) was used. This plate only allows the flow of brine; a high entry capillary pressure prevents the entry of CO 2 . The capillary entry pressure of the plate was checked before the start of the experiment and found to be higher than the pressure drops achieved during the experiment. As the permeability of the porous plate was many orders of magnitude lower than that of the core, practically the entire applied pressure drop occurred across the plate. Before the leading CO 2 -brine interface comes into contact with the porous plate, this pressure drop applies a constant flow rate boundary condition directly at the outlet face of the core, since we only have single-phase (brine) flow through the porous plate. Once the leading CO 2 -brine interface comes into contact with the plate, Detail of the core assembly showing aluminium core wraps to prevent the diffusive exchange of CO 2 across the flouro-elastomer sleeve, a low-permeability porous plate used for flow control, and the PTFE spacer used to monitor the CO 2 -brine interface during system pressurisation. Valves are numbered V 1−14 the brine flow rate through the plate declines. The plate will remain brine-saturated, due to its high capillary entry pressure, and the capillary pressure within the core will increase until it is equal to the applied pressure drop across the plate. Using these boundary conditions removes the impact of the dead volume on the flow through the core, allowing for an accurate flow control at extremely low rates.

Methodology
The brine used contained potassium iodide (KI) with a salinity of 1.5 M. KI was chosen as a solute as it has a high mean atomic number, giving a high X-ray attenuation coefficient. This makes it an effective contrast agent, causing large greyscale differences between the brine and the CO 2 on reconstructed tomographic images. To prevent any reaction during the experiment, either between carbonate rock and carbonic acid or between the brine and the supercritical CO 2 , the three phases were mixed together at working conditions (50 • C, 10 MPa) in a heated reactor (Parr Instruments Co., IL, USA) prior to injection. System pressure was maintained and fluid flow controlled using high pressure syringe pumps (Teledyne ISCO, Lincoln, NE, USA, Model 1000D). Pressure was measured using inbuilt transducers within the pumps, calibrated at high pressures before and after each experiment, with no observed measurement drift.

Experimental Protocol
The drainage experiment consisted of the following steps: 1. The core was imaged along its entire length using a large number of projections (3600) to create a high-quality unsaturated image for both pore-space analysis and comparison with the saturated images. 2. The temperature and pressure in the reactor were raised to that desired for the experiment (50 • C and 10MPa), and the fluids within the reactor were vigorously mixed. 3. The core was then mounted in the coreholder, and the coreholder was brought up to reservoir temperature and pressure (10 MPa, 50 • C), establishing a fully brine-saturated core prior to the experiment. The efficacy of initial saturation of the pore-space could be verified by taking test tomographic scans at this point, and the sample was always found to be 100 % saturated prior to CO 2 injection. During pressurisation, the CO 2 -brine interface can be monitored through the PTFE spacer. 4. CO 2 was injected using a pressure drop across the porous plate of 25 kPa, equating to a single-phase flow rate of 1.75 × 10 −15 m 3 /s. 5. Tomographic images were continually acquired during CO 2 injection, using 800 projections each with an exposure time of 0.04 s. 6. The qualitative progress of the drainage process was monitored without stopping the image acquisition sequence by monitoring changes in the first radiographic projection from each scan. When no further changes were seen in the CO 2 saturation over the course of 30 min, acquisition was stopped.
The image time series was reconstructed using a filtered back-projection algorithm (Titarenko et al. 2010). The reconstruction centre was found for the first and last image in the sequence and linearly interpolated between these two values for the others. Each reconstructed image consisted of around 2200 3 voxels with a voxel size of 1.82 µm, representing a field of view of 4 mm × 4 mm × 4 mm. Imaging was conducted at the I13 beamline of Diamond Light Source under pink beam conditions, with photon energies ranging from 0 to 30 keV. Each image took around 45 s to acquire, with 32 s spent taking projections and around 13 s spent returning to the initial state and preparing for the next scan.

Image Processing
All image processing was conducted using the Avizo Fire 8.0 (Visualization Sciences Group, Burlington, MA, USA) and imageJ programs.

Partially Saturated Images
After reconstruction images were binned, so each 8 (2 × 2×2) voxel cube was averaged to a single voxel to reduce noise and make the data series easier to analyse. The images were then cropped around the core so each one consisted of around 1100×1100×1100 voxels with a voxel size of 3.64 µm. The images were then passed through a 2D non-local means edge preserving filter (Buades et al. 2008(Buades et al. , 2005 and transformed so that any small movement between scans was corrected by the process of registration. This was done by comparing each image to the first in the sequence (the reference image) using a normalised mutual information metric (Pluim et al. 2000(Pluim et al. , 2003 in the process of registration. The registered images were resampled onto the same voxel grid as the first image using a Lanczos filter for voxel-wise interpolation (Duchon 1979). This allowed for multiple transformed images in the displacement sequence to be compared on the same unique grid without concern about later interpolation. The reference image was then subtracted from the resampled image, giving a difference image (Fig. 2c). The use of difference images not only increased the contrast of displaced CO 2 , aiding in segmentation, but also allowed for CO 2 invasion to be imaged using volume renderings without segmentation (Fig. 2d). These volume renderings allowed for the pairs (or sequences) of images associated with each event to be identified for later segmentation and analysis, reducing the computational cost of the imageprocessing workflow. Exemplar events were identified such that all fluid movement occurred between subsequent images, fluid interfaces were well defined, and no fluid movement (as shown by interface blurring or regions with greyscale intermediate between CO 2 and brine) was observed.
Once individual events of interest were identified from the sequence, the data set was segmented into two phases (where CO 2 was treated as one phase and brine and rock were treated together as the single other phase) using the greyscale gradient as an input to a watershed algorithm, seeded using a 2D histogram (Jones et al. 2007). The image-processing workflow is shown in Fig. 2.
Capillary pressure (the pressure difference between the non-wetting and the wetting phase) is related to the curvature of the interface between two immiscible phases by the Young-Laplace equation.
where σ is the interfacial tension,Ĉ is the mean interface curvature, C 1 is the major curvature (=1/R 1 ), and C 2 is the minor curvature (=1/R 2 ) such that C 1 > C 2 (R 2 > R 1 ). Capillary pressure is traditionally measured using techniques such as mercury injection (Purcell 1949), porous plate coreflooding (El-Maghraby and Blunt 2013; Pentland et al. 2011) or quasisteady-state coreflooding (Pini et al. 2012); however, recent developments have allowed for direct measurement of capillary pressure at the pore-scale through the measurement of the curvature of surfaces generated across the wetting-non-wetting phase interface (Andrew et al. 2014a;Armstrong et al. 2012). Surfaces were generated using a smoothed marching cubes algorithm, and curvatures found by creating a best-fit quadratic surface for each of the elements across the CO 2 -brine interface (Hege et al. 1997;Stalling et al. 1998) (Fig. 2g, h) producing a surface scalar field. This method was originally described by Armstrong et al. (2012), and errors in the protocol for the CO 2 -brine-carbonate system addressed by Andrew et al. (2014a). Extending this technique to measure the capillary pressure during drainage in this system is, however, challenging. As drainage capillary pressures are higher than those seen during imbibition-the radii of curvature are smaller (relative to the voxel size)-and the technical measurement of curvature is much more difficult. The CO 2 -brine interface either resides as a corner meniscus, with its curvature pointing into the corners a pore (Fig. 3a, region 2; Fig. 3b, region 2), or resides as a terminal meniscus, with the curvature of the interface pointing into pore-throats (Fig. 3a, region 1; Fig. 3b, region 1). As this system is weakly water-wet (Andrew et al. 2014b), these regions have uniformly positive curvature, making them distinct from the CO 2 -rock interface (Fig. 3a, region 3; Fig. 3b, region 3) which has a negative curvature. At higher capillary pressures, curvature measurements from corner menisci become unreliable as the curvature vectors in pore corners are very anisotropic (R− 1 R 2, Fig. 5e), so although the minor curvature b Images filtered using a non-local means edge preserving filter. c Different images were created by finding the difference between each image and the reference image. This increases the contrast of CO 2 , making segmentation easier and making it possible to view the progress of the drainage front without the need for segmentation of each individual image. d These different images were resampled using an 8 × 8 × 8 kernel, allowing for the entire 4D data set to be examined and rendered in 3D on a single workstation. The region examined in 4G-H is shown as a red cube in 4D. e, f Individual events are identified and the data segmented using watershed seeded using a 2D histogram for the two scans either side of the event. g A surface was generated for both the connected and disconnected CO 2 . h Curvature was found on this surface by fitting a polynomial surface to each of the elements on the CO 2 surface. Each element was then assigned a curvature value from this best-fit surface and an average interface curvature found from resulting distributions from selected areas of the interface surface Curvature distributions across the CO 2 surface can be separated into three different areas, numbered in a-c. Region 1 indicates the portions of the CO 2 -brine interface residing in pore-throats (terminal menisci), region 2 indicates the portions of the CO 2 -brine interface residing in pore corners (corner menisci), and region 3 indicates the CO 2 -grain interface. b Shows a curvature map with the pore-space not occupied by CO 2 shown in red, showing the relationship between terminal menisci, pore bodies and pore-throats. c-d Show the intersection of corner menisci and pore corners. The blue colour represents solid grain. e Shows the curvature vectors resolved on the surface of a corner meniscus and f on a terminal meniscus. The major curvature vector (C 1 ) is shown in red, whereas the minor curvature vector (C 2 ) is shown in black. Corner menisci e show highly anisotropic curvatures (C 1 C 2 ), whereas in terminal menisci they are more isotropic vector (C 2 ) is well resolved (its radius is large relative to the voxel size), the major curvature vector (C 1 ) is very poorly resolved (its radius is smaller relative to the voxel size).
Terminal menisci are, however, better resolved, as they have much more isotropic curvature vectors (R 1 ≈ R 2 , Fig. 5f). The most symmetric curvature vectors (residing at the centres of terminal menisci) are therefore the ones where both R 1 and R 2 are largest relative to the voxel size, meaning the overall measurement of mean curvature (and so capillary pressure) is the most accurate (Fig. 3e, f). The relationship between terminal menisci, pore bodies and pore-throats can be seen visually in Fig. 3b.
To identify these terminal menisci, both curvature vectors (C 1 &C 2 ) were computed from the triangular mesh and a vector isotropy criterion applied such that where F is some contrast factor in the range from 0 to 1. The contrast factor was chosen such that measurements were taken from at least 200 points across the surface, ensuring a representative sampling of the entire interface. A contrast factor of 0.9 (so that measurements were only taken if the magnitude of C 2 associated with any triangular element was 90 % of the magnitude of C 1 for that element) was used for the measurements described in this study.
As the image-processing workflow for each of the images associated with each event was kept constant, resulting changes in curvature must result from changes in fluid distribution and pore-scale capillary pressure. Measurements of the connected CO 2 interface curvature were taken across the entire field of view, apart from Sect. 3.1 where the connected CO 2 phase was split up into 4 equal volume quadrants to test assumptions of capillary equilibrium.
Curvature distributions for the entire exterior CO 2 surface and for terminal menisci (as identified by curvature vector anisotropy, as described above) are shown in Fig. 4. The curvature distribution of the entire CO 2 surface varies continuously from negative to positive, making any unique curvature description associated with the CO 2 -brine interface impossible. The curvature distribution of terminal menisci, however, has a much better resolved uniquely positive curvature. By selecting terminal menisci for measurement using curvature anisotropy, we can effectively maximise the range of capillary pressures accessible to measurement by curvature analysis, allowing us to measure higher capillary pressures. What is more, as the regions of the interface selected by the technique are objective, any changes in the distribution are likely to be caused by changes in capillary pressure. Also, these measurements can be made from two-phase segmentation, rather than (much more difficult) three-phase segmentation.
To find an objective estimation of peak position and width, these distributions were then fitted to a Gaussian model using a trust region algorithm (Conn et al. 2000). This peak position was then converted to a pressure difference across the interface using the Young-Laplace equation (Eq. 1). The interfacial tension in this system is estimated by linearly interpolating between measurements for a given pressure, temperature and salinity (Li et al. 2012;Georgiadis et al. 2010), and a value of 0.037 N/m was used in this study.
As measurement biases can never be completely eliminated, results are reported not only using histograms of the curvature of terminal menisci, but also using curvature maps over the entire ganglion within the field of view (e.g. Fig. 4a, b) to allow for visual verification of the resulting conclusions. From these curvature maps, and from the resulting measured curvature distributions, we can draw conclusions about flow phenomena based on the relative changes in curvature for connected and disconnected regions of CO 2 .
The fitted curvature for the ganglion shown in Fig. 4 is 0.027±0.004 µm −1 , corresponding to a radius of curvature of 37 ± 6 µm (around 10 times larger than the voxel size). Fig. 4 Selective measurement of mean interface curvature. a Shows the interface curvature distribution for an entire connected CO 2 ganglion. b Shows the terminal menisci, as selected using curvature anisotropy. These regions were used for curvature measurement. c Shows a comparison of the curvature distributions for measurements taken from the entire ganglion (blue) and for measurements taken only from the terminal menisci (red). Curvature measurements taken only from the terminal menisci show a single well-defined peak, giving the most accurate representation of CO 2 -brine interface curvature

Dry Scans
The dry scan was reconstructed, binned and filtered in the same way as the partially saturated images. The dry scan was then registered to the saturated reference scan using the same registration algorithm described above. Finally, the reconstructed and registered dry scan was resampled onto the same grid as the partially saturated images using a Lanczos filter (Duchon 1979). The dry scan was segmented into pore and solid using a watershed algorithm seeded using a 2D histogram. The pore-space was separated into individual pores and throats by first calculating a Euclidean distance map of the pore-space. A separate pore was then defined for each of the catchment basins generated by running the watershed algorithm on this distance map (Beucher and Lanteujoul 1979;Wildenschild and Sheppard 2013). Throats were created across the interfaces between different catchment basins. Pore and throat radii were defined as the maximum value of the Euclidean distance map function for each pore or throat, equivalent to the radius of a maximally inscribed sphere within the element (Fig. 5).

Results
Drainage does not progress as a uniform front, but instead is quantised into a discrete set of events where regions of the pore-space occupied by brine in one image are filled with CO 2 in the next. The interface movements due to these jumps are not resolved with the temporal resolution available for this method (45 s); micromodel studies have shown that pores drain on millisecond timescale regardless of the (externally imposed) capillary number (Armstrong and Berg 2013;Mohanty et al. 1987). This invasion mechanism supports the conclusion that, in this system, CO 2 is the non-wetting phase, as shown in Andrew et al. (2014a). Events imaged using this technique represent changes between two apparent equilibrium states, before and after each of the non-wetting phase jumps. If a jump occurs during the acquisition period of a scan, then the reconstructed greyscale of the volume corresponding to the jump represents a time-averaged value between the two fluid phases.
A typical jump can be seen in Fig. 7. Between two scans, a region of the pore-space previously occupied by brine becomes occupied by CO 2 , connected with the rest of the connected CO 2 phase through the throat invaded during the jump. As the time resolution of this study (45 s) is limited, it is impossible to say whether these events are the result of a single large Haines jump or a sequence of rapidly occurring smaller Haines jumps. However, it is still possible to draw conclusions on flow phenomena based on relative changes in fluid position and interface curvature.

Capillary Equilibrium
First, we will examine the assumption that between rapid drainage events the interface is at capillary equilibrium. To do this, the connected CO 2 is split into four equal volume quadrants, and the curvature is calculated across the entire surface within each quadrant to identify the terminal menisci using curvature anisotropy. The resulting distribution of curvature in the terminal menisci for each quadrant, as well as when summed across the entire interface is shown in Fig. 6. These measurements were taken from an image that shows no evidence of a rapid drainage event during the scan acquisition.
The resulting curvature distributions from each quadrant coincide, supporting the conclusion that, at the spatial and temporal scales examined with this technique (mm and tens of seconds, respectively), there are no gradients in capillary pressure across the interface. This, combined with the observation that CO 2 fills regions of the pore-space rapidly (between one scan and the next), supports the findings of Mohanty et al. (1987) and Armstrong and Berg (2013) that fluid displacement occurs at much shorter timescales (at most tens of ms).
Displacement sequences occurring during the scan time would be represented as blurred regions in the reconstructed data set acquired during the time that the jump occurred. The Image-processing workflow for dry scans of the pore-space. a The image after reconstruction and filtering. b The image was segmented using a 2D histogram-based watershed segmentation (Fig. 2). c A Euclidean distance map of the pore-space was calculated by assigning each voxel the distance from it to the nearest pore wall. d This was then used to separate the pore-space out into individual pores by calculating the watershed basins of this distance map. Pore-throats are defined as the local minima in the Euclidean distance map, corresponding to constrictions in the pore-space-they are surfaces dividing pores. The centres of pore bodies represent the local maxima of the distance map. e-f This separation process is shown for two example pores. e Regions of high Euclidean distances (shown in red) are defined as pore bodies and throats are defined as the surfaces connecting them. f The two separated pores from the distance map shown in e drainage events examined below in Sects. 3.2-3.4 were chosen not to show this effect, indicating that displacements occurred between the end of one scan and the start of the next.

Capillary Pressure Change and Interface Recession
The first drainage event examined in this study is shown in Fig. 7.
The jump causes the interface residing in the corners of pores and pore-throats to recede associated with a change in capillary pressure (as shown by the red regions in Fig. 7). This  The region occupied by CO 2 after the jump but not before is shown in blue; the region occupied by CO 2 in both scans is shown in yellow; the region occupied by CO 2 before the jump but not after is shown in red. The quasi-static capillary pressure change associated with the event causes the CO 2 -brine interface to recede in pore-throats and the corners of pores. b The same rendering of a subvolume around the event, with the subvolume taken from the region indicated by the blue cuboid in a recession provides the fluid required for the jump, allowing for rapid pore filling without a discontinuity in bulk fluid volume. Changes in capillarity associated with drainage events are, therefore, non-local, contrary to the quasi-static assumptions present within many network models. Instead, each jump affects fluid configurations across the entire wetting-non-wetting phase interface. The resulting curvature maps can be seen in Fig. 8, and a schematic representation of this process is shown in Fig. 9.
In the example displacement studied, the fitted curvature peak position prior to the jump is 0.033 ± 0.003 µm −1 corresponding to a capillary pressure of 2.45 ± 0.2 kPa and after the jump is 0.0268±0.0015 µm −1 , corresponding to a capillary pressure of 1.98±0.11 kPa. The non-wetting phase residing in pore-throats across the CO 2 -brine interface supplies the fluid  Mohanty et al. (1987) and Yuan (1990) and micromodel studies by Armstrong and Berg (2013), providing the first evidence for this process by direct pore-scale observation in a rock at reservoir conditions.

Roof Snap-Off
Drainage events do not always proceed as a single drainage event. As the wetting-nonwetting phase interface emerges from a pore-throat during a drainage event, the curvature across the leading interface of the non-wetting phase will decrease. This will lead to curvature disequilibrium between different regions of the interface. As the interface attempts to reequilibrate, the curvature in the pore-throat connecting the leading interface to the rest of the non-wetting phase will decrease. If it decreases below the snap-off capillary pressure of Fig. 9 A schematic representation of a drainage event. The non-wetting phase (CO 2 ) is shown in red, the schematic solid grains are white, and the wetting phase (brine) is shown in blue. When the interface capillary pressure overcomes the threshold invasion capillary pressure for throat 1, the next pore is rapidly invaded, and a new, more stable state established with the next pore saturated with the non-wetting phase. The capillary pressure of the new state will be lower than prior to the jump, causing the wetting-non-wetting phase interface to recede in the throats labelled 2 that throat, the leading portion of the non-wetting phase will disconnect, forming an isolated ganglion in the process of snap-off.
Although snap-off is more commonly associated with macroscopic imbibition (wetting phase ingression), it was first described experimentally and analytically for the process of drainage by Roof (1970), with further recent quantitative analysis using finite volume modelling by Raeini et al. (2014a). The mechanism for this process can be seen schematically in Fig. 10.
Roof snap-off occurs when the leading interface advances through the second pore, reducing its interfacial curvature. The surface will try to maintain capillary equilibrium, causing the wetting layer in the throat (Fig. 10, throat 1) to thicken until it becomes unstable, at which point wetting phase will rapidly fill the throat, disconnecting a droplet of the advancing non-wetting phase from the rest of the connected non-wetting phase. After this the disconnected ganglion may rearrange, as it moves towards the centre of the pore in which it resides, potentially changing its capillary pressure.
This phenomenon can be seen in an example sequence of four successive images within our data set, each separated by 45 seconds (Fig. 11). Surface curvature was measured, as described in Sect. 2, for both the connected and disconnected CO 2 phase, and the results can be seen in Fig. 11.
Between steps 1 and 2 (Fig. 11a, b), a small disconnected droplet (ganglion) of CO 2 with a volume of 4.33 × 10 6 µm 3 appears in the pore-space, due to a cycle of local invasion and snap-off. Between steps 2 and 3 (Fig. 11b, c), another cycle-non-wetting phase invasion and a second snap-off event-causes this ganglion to grow to a volume of 1.18 × 10 7 µm 3 [multiple cycles of this process are not unexpected and were shown in direct simulation studies by Raeini et al. (2014a)]. Finally, by step 4 it is connected with the remainder of the CO 2 phase. The curvature distributions for the connected and disconnected ganglia are shown in Fig. 11e. A summary of the fitted curvature peak positions and corresponding capillary pressures can be seen in Table 2.
From the curvature distributions of the connected and disconnected non-wetting phase, we can observe two interesting features. Fig. 10 A schematic of the process involved in Roof snap-off, where a ganglion becomes disconnected immediately after throat invasion by snap-off of the throat through which the invasion took place. a Shows the situation before the jump. During the event, extremely low interface curvatures are generated locally around the event, causing the wetting layer in throat 1 to swell (b), eventually leading to its snap-off and the disconnection of non-wetting phase in an isolated ganglion (c). After the ganglion has become disconnected, there may be some subsequent ganglion rearrangement as a local surface energy minimum is found First, the curvatures of the disconnected ganglion after the first and second cycle of local invasion and snap-off are very similar despite the ganglion almost trebling in volume after the second snap-off cycle. Their similarity indicates that the primary control on the capillary pressure of residual ganglia is not the rearrangement of the disconnected ganglion after snap-off, but the threshold snap-off capillary pressure of the disconnecting throat. This is independent of ganglion volume as it is defined by the local topography of the pore-throat. We would expect the impact of ganglion rearrangement on capillary pressure to be different for ganglia of different volumes, with small ganglia rearranging more than large ganglia. This may be the cause of the slight difference in curvature distribution between the first and second cycle of Roof snap-off [with the capillary pressure of the ganglion present in step 3 (Fig. 11c, e curve v) being slightly higher than that present in step 2 (Fig. 11b, e curve vi)]. This change is close to the measurement error, much smaller the difference between the curvature of the connected and disconnected fluid interface. This finding supports the findings of Andrew et al. (2014a) correlating the capillary pressures of a set of ganglia at residual saturation with the inverse of local throat radius (and so snap-off capillary pressure), a correlation which should only be present if the snap-off capillary pressure of the disconnecting throat defines the capillary pressure of the disconnected ganglion.
Second, the capillary pressures of the disconnected ganglia at both time steps are much lower than the capillary pressure of the connected CO 2 at any point during the sequence. We propose that this is due to significant dynamic capillary pressure gradients generated in the leading region of the non-wetting phase during the drainage event, out of capillary equilibrium. These cause the low capillary pressures required for throat snap-off, which are then preserved in the disconnected ganglion. Snap-off is not caused, therefore, by the Fig. 11 a-d Curvature maps for the connected and disconnected CO 2 regions, rendered in a subvolume from four sequential images during a Roof snap-off sequence. Although a subvolume is used in a-d for clarity, curvature measurements for the connected CO 2 were made across the entire CO 2 surface. e The resulting capillary pressure distributions. Low capillary pressures in the disconnected phase region may be result of low transient capillary pressures generated during the drainage event Table 2 Summary of fitted curvature peaks and corresponding capillary pressures for both the connected and disconnected wetting phase during the event seen in Fig. 11 Step Connected curvature (µm −1 )  Fig. 12 A rendering of a large drainage event, coloured using the same scheme as Fig. 6, with the regions occupied by CO 2 after the jump but not before shown in blue, the regions occupied by CO 2 both before and after the event shown in yellow and the regions occupied before but not after the event shown in red. This shows not only the interface recession resulting from the quasi-static capillary pressure change associated with the drainage event, but also the snap-off of one of the throats occupied by CO 2 prior to the jump. Unlike the event analysed in Sect. 3.3, this throat is not that through which the invasion occurred during the drainage event. Figure 13 shows local renderings of this snap-off event with the capillary pressure evolution for both the connected and the disconnected phase longer timescale capillary pressure changes occurring across the entire surface, as examined in Sect. 3.2. After the event, the rest of the connected CO 2 -brine interface then rapidly equilibrates to capillary pressures close to those seen prior to the jump (on time scales shorter than those accessible using this imaging technique).

Distal Snap-Off
Snap-off associated with drainage events is observed not only locally, in the throat invaded by the event, but also further away from the event location. This non-local process is termed "distal snap-off", and an example is shown in Fig. 12. The curvature of the connected and disconnected non-wetting phase regions before and after the drainage event is shown in Table 3.
As the disconnected ganglion curvature is lower than the connected CO 2 curvature both before and after the drainage event, the observed snap-off cannot be explained by the measured Fig. 13 a, b Non-wetting phase interface curvature maps of the volume surrounding the snapped-off throat shown in Fig. 12 before the jump (a) and after it (b). Although a subvolume of the curvature map is shown for clarity, connected CO 2 curvature measurements were taken from the entire CO 2 surface. c Shows the curvature distributions of the terminal menisci, indicating that the ganglion curvature (and so capillary pressure) is much lower than the curvature (and so capillary pressure) of the connected CO 2 both before and after the event  .2). Instead, the snap-off must be caused by some transient low dynamic capillary pressure seen during interface movement, which is then preserved in the capillary pressure of the disconnected ganglion. During a large jump, such as that seen here, the non-wetting and wetting phase are displaced and consequentially interface curvature is reduced, with these changes propagating away from the leading interface. This causes not only the interface to recede from the throats and corners directly involved in the non-wetting phase invasion, but can potentially reduce the capillary pressure in throats a significant distance away. This can reduce the capillary pressure in these throats to below their snap-off capillary pressure, disconnecting regions of the non-wetting phase (Fig. 14).
The controls on snap-off in arbitrarily complex geological media are difficult to understand. However, it is useful to consider the case of the hypothetical geometric elements (typically with triangular cross sections) present in network models (Valvatne and Blunt 2004). As capillary pressure decreases in these systems, the wetting phase present in the corner regions of elements swells until they come into contact. At this point, the fluid-fluid interface becomes unstable as the curvature required to maintain the interface is higher than those in the surrounding elements, causing wetting phase to rapidly fill the element in the process of snap-off. If the wetting phase in all three corners of a triangular element is swelling concurrently (relying on the assumption that the wetting phase is well connected) snap-off will occur at a capillary pressure of: where θ a is the advancing contact angle and β is the effective throat corner half-angle (Valvatne and Blunt 2004). The snap-off capillary pressure is therefore a function of throat size and shape, with smaller throats snapping off at higher capillary pressures than larger ones, and throat radius being the primary control on pore-throat snap-off capillary pressure (Andrew et al. 2014a). The rest of the throats occupied by the non-wetting phase at the time of the jump remain occupied after the event has occurred. Although modelling the flow through each of the throats is computationally challenging, as the primary control on snap-off capillary pressure is throat radius (Sect. 3.3), we can use the throat radii (determined as described in section 2.2.2.2) to examine the behaviour of other throats within the connected CO 2 ganglion.
An example throat was identified which remained saturated with CO 2 throughout the event, despite having a smaller radius than that of the snapped-off throat (snap-off occurred in a throat of radius 26 µm, whereas a throat with a radius of 19 µm remained saturated with CO 2 ). As throat radius is the primary control on snap-off capillary pressure, and both were subject to the same dynamic forces during the drainage event, we might expect for snap-off to occur in both of these throats, or at least in the smaller one, rather than the other way round.
One possible explanation is proximity. We would expect the magnitude of dynamic changes in capillary pressure to decrease further away from any event; however, the narrow throat which remained saturated with CO 2 was much closer to the event than that which snapped off (440 µm as opposed to 1640 µm), so it would be subjected to similar or larger dynamic forces than those experienced by the snapped-off throat.
One critical difference may be the availability of wetting phase and how quickly it can invade the throat. The smaller throat which remains full of non-wetting phase after the event does not connect to any pores which are adjacent to pores saturated with the wetting phase. The only path through which the wetting phase can invade the throat is therefore through wetting layers along roughness or corners in the solid surface. Although these may be present, the rate at which the wetting phase can move through these layers to a throat deep within the interior of the connected non-wetting phase may be insufficiently rapid to allow snap-  Fig. 12. A Haines jump of fluid into pore 1 not only causes interface recession, but the low dynamic pressure generated during the event is sufficient to cause snap-off in the throat 3. Although throats 2 and 3 have very similar snap-off capillary pressures, throat 2 does not snap-off as it is further away from available wetting phase, meaning that the wetting layers could not grow sufficiently for snap-off to occur in the fast dynamic timescales involved in the jump process off on the ms timescale associated with the jump. As the longer timescale capillary pressure change (as examined in Sect. 3.2) is not sufficient to cause snap-off these throats, they remain saturated with CO 2 .
In contrast, the throat that did snap-off is connected to pores which were adjacent to pores filled with the wetting phase, allowing more wetting phase flow during the event. This then allowed dynamic pressure gradients created during the event to be reflected by wetting layer growth within the throat, causing snap-off (Fig. 14, throat 3).
The process of distal snap-off is potentially more important for the long-term evolution of the drainage sequence than local Roof snap-off (Sect. 3.3) as the fluid configurations arising from distal snap-off are more persistent. The throat through which any specific Haines jump occurs is selected because it is the most energetically favourable; it has the lowest entry capillary pressure for any of those accessible to the wetting-non-wetting phase interface. If snap-off occurs through the throat through which the drainage event occurs without significant rearrangement in the connected non-wetting phase, then this throat will remain the most energetically favourable for a filling event. It will therefore be the site of the next event in the sequence. The disconnected ganglion volume will gradually increase due to a sequence of local invasion snap-off events until it becomes connected with the rest of the (connected) non-wetting phase. If, however, snap-off occurs in another throat unrelated to the throat through which the drainage event occurs, there is no reason that the snapped-off throat should be refilled until the capillary pressure increases to above its invasion threshold. The disconnected ganglion can therefore remain disconnected for much longer in the drainage sequence, potentially impacting the wetting phase flow field and further dynamic events. This is seen by comparing the jump discussed in this section to that discussed in Sect. 3.3. Whereas the disconnected ganglion formed due to the jump discussed in Sect. 3.3 is reconnected fairly quickly with the rest of the connected CO 2 (within the time needed to acquire 4 tomographic images, or around 3 min), the disconnected ganglion seen in Fig. 13 is not reconnected with the rest of the non-wetting phase for the time needed to acquire another 73 images (around 55 min). This is much later in the drainage sequence, where the leading non-wetting phase interface had moved out of the field of view, meaning that this ganglion (preserving features associated with the dynamics of non-wetting phase jumps) has a persistent impact on the wetting phase flow field and the sequence of subsequent drainage events.

Conclusions
We have used synchrotron-based fast X-ray microtomography to image drainage in a carbonate sample at representative subsurface conditions with a time resolution of around 45 s. Capillary pressure evolution around these events has been calculated using curvature mapping of terminal menisci, identified using curvature anisotropy. By focusing on a small number of discrete drainage events in detail and comparing capillary pressure measurements to localised multiphase fluid modelling, we have observed dynamic drainage phenomena at representative reservoir conditions.
Quasi-static capillary pressure changes associated with these events were found and measured. We provide the first evidence for non-local interface recession associated with a change in capillary pressure in a real rock at reservoir conditions. Furthermore, snap-off was observed, both in the throat through which non-wetting phase invasion occurred (Roof snap-off) and in throats far away from it (distal snap-off). It was found that distal (non-local) snap-off produced more persistent fluid configurations than Roof snap-off. Equilibrium capillary pressure changes associated with drainage events are not sufficient to explain the snap-off associated with these events, as the disconnected ganglia formed have lower capillary pressures than the connected non-wetting phase either before or after the jump. The snap-off may instead be due to significant dynamic pressure gradients generated during these events as fluid pressures and interface curvatures attempt to re-equilibrate after rapid fluid movement. Snap-off in these throats is not only controlled by throat radius and proximity to the event but also by the local fluid arrangement. A more mobile local wetting phase may lead to a greater likelihood for a throat to snap-off during an event.
Future work could focus on examining the relative statistical importance of each of the phenomena observed and identified in this paper to the overall drainage process. Another interesting avenue of future research would be the interaction between macroscopic flow rate and the length scale over which non-local capillarity is observed in larger systems. The method described in this paper for controlling extremely low flow rates at subsurface conditions may be of use in future work using dynamic tomography to address many other issues in multiphase flow in porous media.