Accepted Manuscript

Connected Pathway Relative Permeability from Pore-Scale Imaging of Imbibition

PII: DOI:

Reference:

S. Berg , M. Rucker, H. Ott, A. Georgiadis , H. van der Linde , F. Enzmann , M. Kersten , R.T. Armstrong , S. de With , J. Becker, A. Wiegmann

S0309-1708(16)30011-2 10.1016/j.advwatres.2016.01.010 ADWR 2549

To appear in:

Advances in Water Resources

Received date: 25 June 2015 Revised date: 25 January 2016

Accepted date: 30 January 2016

Please cite this article as: S. Berg , M. Rucker, H.Ott, A. Georgiadis , H. van der Linde , F. Enzmann , M. Kersten , R.T. Armstrong , S. de With , J. Becker, A. Wiegmann , Connected Pathway Relative Permeability from Pore-Scale Imaging of Imbibition , Advances in Water Resources (2016), doi: 10.1016/j.advwatres.2016.01.010

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

Highlights

• Quasi-static simulation of two-phase flow in porous media agrees with experimental data within experimental uncertainty for drainage

• A morphological approach, which approximates capillary displacement, does not represent the imbibition process. Ultimately for modelling relative permeability in imbibition an approach is needed that captures moving liquid-liquid interfaces which requires viscous and capillary forces simultaneously

• If pore scale fluid distributions are available e.g. from micro-CT flow experiment relative permeability can be estimated from the connected pathway flow (for low capill / numbers)

• The agreement is better at low water saturations where the oil phase is predom antly connected than at higher water saturation where the oil phase is increasingly di Dnnected.

Connected Pathway Relative Permeability from Pore-Scale Imaging of Imbibition

S. Berg1-*, M. Rücker1'2, H. Ott3,1, A. Georgiadis1, H. van der Linde1, F. Enzmann2, M. Kersten2, R. T. Armstrong4, S. de With5, J. Becker6, A. Wiegmann6

1 Shell Global Solutions International B.V., Kessler Park 1, 2288 GS Rijswijk, The Netherlands

2 Geosciences Institute, Johannes Gutenberg University, Becherweg 21, 55099 Mainz, Germany 3Department of Petroleum Engineering, Montanuniversität Leoben, A-8700 Leoben, Austria 4School of Petroleum Engineering, University of New South Wales, NSW 2052 Sydney, Australia 5Technical University of Delft, The Netherlands 6Math2Market GmbH, Stiftsplatz 5, 67655 Kaiserslautern, Ger

"corresponding author: steffen.berg@shell.com

Abstract

Pore-scale images obtained from a synchrotron-based X-ray computed micro-tomography (^CT) imbibition experiment in sandstone rock were used to conduct Navier-Stokes flow simulations on the connected pathways of water and oil phases. The resulting relative permeability was compared with steady-state Darcy-scale imbibition experiments on 5 cm large twin samples from the same outcrop sandstone material. While the relative permeability curves display a large degree of similarity, the endpoint saturations for the ^CT data are 10% in saturation units higher than the experimental data. However, the two datasets match well when normalizing to the mobile saturation range. The agreement is particularly good at low water saturations, where the oil is predominantly connected. Apart from different saturation endpoints, in this particular experiment where connected pathway flow dominates, the discrepancies between pore-scale connected pathway flow simulations and Darcy-scale steady-state data are minor overall and have very little impact on fractional flow. The results also indicate that if the pore-scale fluid distributions are available and the amount of disconnected non-wetting phase is low, quasi-static flow simulations may be sufficient to compute relative permeability. When pore-scale fluid distributions are not available, fluid distributions can be obtained from a morphological approach, which approximates capillary-dominated displacement. The relative permeability obtained from the morphological approach compare well to drainage steady state whereas major discrepancies to the imbibition steady-state experimental data are observed. The morphological approach does not represent the imbibition process very well and experimental data for the spatial arrangement of the phases are required. Presumably for modelling imbibition relative permeability an approach is needed that captures moving liquid-liquid interfaces, which requires viscous and capillary forces simultaneously.

1 Introduction

For the description of multiphase flow in porous media on the Darcy scale, commonly the (phenomenological) extension of Darcy's law to multiphase is employed. The concept of relative permeability is introduced (Dullien, 1979; Corey, 1954) to account for the presence of other immiscible fluid phases (Niessner, 2011). Together with the standard parameters porosity (absolute) permeability K, fluid viscosities u, and the capillary pressure pc vs saturation 5W function, the relative permeability kr is one of the key parameters to model the flow of multiple phases. For instance, in reservoir engineering, where on the field scale viscous forces dominate over capillary forces (Hilfer, 1996), the relative permeability vs saturation function together with absolute permeability has a large influence on the flux of fluid phases. Relative permeability can largely vary depending on rock type, wettability and other parameters, which can have a major impact on oil recovery by water flooding. At the same time, relative permeability is difficult and is experimentally time-consuming to obtain.

Since the multiphase extension of Darcy's law is phenomenological it does not within its own boundaries predict relative permeability. Relative permeability is traditionally determined experimentally using either the steady state (Dake, 1978; Kokkedee, 1996) or non-steady state methods (Johnson, 1959; Berg, 2013b), which are time-consuming and costly. Therefore within the last decade there has been increasing interest in predicting relative permeability by pore-scale simulation (Blunt, 2002), e.g. using pore network modelling or direct simulation techniques (Ferrari, 2013, Koroteev, 2014; Landry, 2014; Gray, 2015). Generally these approaches can be categorized into pore network models (Lenormand, 1988 Raoof, 2013; Ghanbarzadeh, 2014) and direct approaches depending on whether the modelling is performed on a network model abstracted from the pore space or directly on the pore space obtained from X-ray computed micro-tomography (uCT) (Prodanovic, 2006; Porter, 2009; Landry, 2011; Wildenschild, 2013; Landry, 2014; Manuel, 2014; Fusseis, 2014).

Pore-scale modelling approaches can be further categorized into quasi-static and dynamic approaches where capillary and viscous forces act at the same time (Leverett, 1941). Combinations of all categories exist, i.e. there are quasi-static (Blunt, 2002) and dynamic pore network models (Joekar-Niasar, 2012; Hammond, 2012; Prodanovic, 2014), but also quasi-static (Becker, 2008) and dynamic direct simulation methods (Coles, 1998; Jaqmin, 1999; Raeini, 2014a,b; Ferrari, 2013; Ramstad, 2009, Vogel, 2005; Demianov, 2011; Boek, 2014; Gray, 2015). In practice, the most common approaches are quasi-static pore network models, which are computationally less costly than dynamic and direct approaches due to simplifications made. However, these simplifications introduce uncertainty. First, the abstracted pore networks are often subject to tuning on reference data (e.g. capillary-pressure data), which decreases a direct link between the actual pore space and the abstracted network (Raeesi, 2013). Elementary pore-scale displacement events (Lenormand, 1988) are captured only in a mechanistic way. Then the capillary-viscous flow is decomposed into capillary-dominated steps determining the pore-scale saturation distribution and viscous-dominated steps where the hydraulic conductivity of connected phases are determined. The main limitation of quasi-static pore network modelling is that only the connected pathway flow is captured. However, multiphase flow in porous media consists of connected flow and ganglion dynamics (Youssef, 2014; Rucker, 2015) caused by the motion of disconnected fluid phases (Mason, 1983; Yadav, 1983; Datta, 2014). In reality there is mass exchange between connected and disconnected fluid phases (Hilfer, 1998; Hilfer, 2006a; Hilfer, 2006b), which not only adds to the net flux but also influences the time evolution of the pore-scale fluid distributions, including the configuration of connected phases.

Dynamic pore network models are more sophisticated and attempt to capture dynamics effects like capillary-viscous displacement including ganglion dynamics correctly (Joekar-Niasar, 2012; Hammond, 2012), but are more expensive from a computational perspective. Direct dynamic simulation approaches based on the lattice Boltzmann technique (Coles, 1998; Ramstad, 2009;

Chukwudozie, 2013; Boek, 2014; Gray, 2015) and Navier-Stokes approaches (Jacqmin, 1998; Demianov, 2011; Raeini, 2014a, Raeini, 2014b; Ferrari, 2013) are becoming more popular because no a priori choice between level of rigor and captured phenomena are made. Like other dynamic approaches, capillary and viscous forces act at the same time and, depending on the choice of flow parameters, both capillary and viscous dominated flow are captured (Raeini, 2014b). That allows for the description of a wide range of flow regimes and simulation of a wide range of pore-scale dynamics such as cooperative and/or non-local displacement processes (Holtzman, 2015). Such nonlocal processes have been observed during drainage in 2D micromodels (Armstrong, 2013; Armstrong, 2015) and during imbibition (Berg, 2015) in with fast ^.CT (Berg, 2013a; Berg, 2013b; Andrew, 2015; Youssef, 2014). However, from a computational perspective they are much more expensive than pore network models.

Along with the choice of the size of the computational domain and associated length scale, often implicitly a choice on the underlying flow regime is also made. Avraam and Payatakes (Avraam, 1995) demonstrated in 2D micromodels during fractional flow that a ganglion dynamics flow regime exists in combination with connected pathway flow. For ganglion dynamics the non-wetting phase is transported in the form of non-percolating clusters (Youssef, 2014). These move through the pore space in a sequence of breakup and coalescence processes (Rücker, 2015), which requires direct numerical simulation to capture (Demianov, 2011; Gray, 2015). Clusters are ultimately mobilized directly or indirectly by viscous forces. Direct mobilization occurs e.g. during "capillary de-saturation" (Mason, 1983; Yadav, 1983; Lake, 1989) where the viscous pressure drop caused by the average wetting phase flow field over non-wetting phase clusters exceeds the capillary trapping forces. For field-relevant flow rates such non-wetting phase clusters can range from a few micrometers (Iglauer,

2012) to centimeters (Armstrong, 2014), which easily exceeds the simulation domain size of direct dynamic simulations. Therefore in a respective numerical experiment aiming to determine relative permeability, the contribution of large non-wetting phase clusters to the overall non-wetting phase flux is commonly dismissed.

Mobilization of non-wetting phase clusters can occur in a more subtle way, which has been discovered in a recent study (Rücker, 2015). Cluster dynamics can also occur in capillary-dominated flow, i.e. where the viscous forces associated with the average flow field are not sufficiently large enough to overcome capillary forces and thus cause mobilization. However, the local flow field velocities associated with pore-scale displacement events can exceed the average flow field velocities by orders of magnitude. Haines jumps and snap-off processes are relatively fast events (Armstrong,

2013) reaching Reynolds numbers near unity (Armstrong, 2015). These events range over many pores, indicating that besides capillary and viscous forces also inertial forces become relevant (Ferrari, 2013) to describe the pore-scale displacement processes. In this way, snap-off processes can trigger coalescence in distant pores connected to the same non-wetting phase cluster (Rücker, 2015), which adds to the viscous flow through connected pathways for the non-wetting phase flux. These processes occur on much shorter length scales than the one defining the balance of viscous and capillary forces (Armstrong, 2014). The associated processes are typically sufficiently captured by direct dynamic simulators using relatively small computational domain sizes (Koroteev, 2014; Raeini, 2014a,b).

Our hypothesis is that if all of the relevant length and time scales associated with the possible flow regimes are captured correctly, then a 1-to-1 match between experiment and simulation would become possible. However, for the computation of relative permeability the question arises of whether there is a significant difference in practice. Large clusters mobilized by viscous forces do contribute to the non-wetting phase flux, no doubt. But if they are infrequent (as percolation models suggest, see Iglauer, 2011 and Georgiadis, 2013), then their contribution may become negligible. The same holds true for the cluster dynamics caused by non-local and cooperative displacement on smaller length scales. To answer this question, we compared Darcy scale steady-state relative permeability to the relative permeability computed (Mostaghimi, 2013) from pore-scale fluid

distributions obtained from a synchrotron-based ¡CT imbibition experiment (Berg, 2015; Rucker, 2015). This approach was already suggested by (Turner, 2004; Hussain, 2014), who performed quasi-static drainage ¡CT experiments, where flow was stopped for imaging, and computed relative permeability from the connected pathway fraction of the resulting fluid distributions using a Stokes flow solver. We advanced this idea in several ways. First, we used pore-scale fluid distributions from a dynamic synchrotron-based ¡CT imbibition experiment (Berg, 2015; Rucker, 2015). Dynamic means that fluid distributions were imaged under flowing conditions, i.e. the visco-capillary balance was maintained. Second, while in previous studies the approach was used for drainage processes, we performed an imbibition study on sandstone rock. From many pore network modelling studies it is known that imbibition processes are less approachable by quasi-static modelling approaches than drainage processes (Mason, 2013). While quasi-static approaches consider only the connected pathway contribution, in reality snap-off (Raoof, 1970) and other mechanisms cause the non-wetting phase to be at least partially disconnected. Therefore, connected pathway flow and ganglion dynamics may co-exist (Rucker, 2015). For validation of the pore-scale modelling, we then compared steady-state relative permeability measurements conducted on a larger-dimension twin sample from the same block to the pore-scale simulation data from the ¡CT experiment. Lastly, to assess whether quasi-static modelling is also sufficient to describe the pore-scale fluid distribution, i.e. for situations where pore-scale fluid distributions are a priori not available, a morphological approach (Becker, 2008; Hilpert and Miller, 2001; Vogel, 2005; Silin, 2010) was used to simulate the capillary-dominated invasion process.

2 Materials and Methods

2.1 Porous Rock

Gildehauser sandstone (Dubelaar, 2015) was used, which has an average porosity of 20% and a permeability of 1.5 ± 0.3 D. It has low clay content and a narrow pore-size distribution (<35 um pore throat diameter), which makes the rock well suitable for ¡CT imaging . Also, the relatively high permeability and low capillary entry pressure make the rock suitable for saturation with a low-pressure flow apparatus.

2.2 Pore-5cale Flow experiments and Pore-5cale Imaging

Cylindrical Gildehauser rock samples of 4 mm diameter (which is larger than the representative elementary volume - REV, see Figure A1 in Appendix A) and 10 mm length were embedded into a cylindrical polycarbonate flow cell by heat-shrinking. The assembly was mounted on a micro-pump assembly, which includes a micro-piston pump for fluid injection. More details about the setup are described in (Berg, 2013a). N-decane was used as non-wetting phase and brine doped with cesium chloride (10 wt-%) as the wetting phase. The interfacial tension between both fluids was a = 35 mN m-1.

The rock sample was first saturated with brine and then n-decane was injected at high rate (1000 ^L/min) to reach a non-wetting phase saturation around 78%. Then the imbibition experiment was initiated by injecting brine at a volumetric injection rate of 0.1 ^L/min, corresponding to a

microscopic capillary number of Camicro = 1.8-10-8 which is defined as Ca

= VWV/a (Lake, 1989).

Here ¡/=1 mPa-s is the wetting phase viscosity, v the linear flow velocity of the injected wetting

macro _ 7 c/ / p i

phase. The cluster-based or macroscopic capillary number Ca = / / r CUmicro (Armstrong, 2014; Berg, 2015) was between 9.1-10-7 and 1.3-10-6 (lcl is the length of the non-wetting phase clusters determined from the ¡CT images and f is the average pore throat radius), indicating that the flow regime (based on average flow rates) is capillary-dominated. During the flow experiment, 60% of the non-wetting phase (in total 2.5-109 ^m3) was mobilized. The flow experiment is described in more

detail in (Berg, 2015; Rücker, 2015), where saturation as a function of time and the (homogeneous) saturation distribution(s) are also shown.

During imbibition the pore-scale saturation distribution was imaged by synchrotron-based fast X-ray micro-tomography under dynamic flow conditions at the TOMCAT beamline at the Swiss Light Source, Paul Scherrer Institute, Villigen, Switzerland. The beam energy was 37 keV (monochromatic, with a multilayer filter). For flow experiments, 1500 projections were collected at 20 ms acquisition time using an edge camera (PCO AG, Kelheim, Germany). The spatial resolution was 2.2 ^m/voxel and the temporal resolution for a 3D rendering was 45 s. Scans of the dry samples were obtained at 60 ms acquisition time to obtain a higher-quality image, as displayed in Figure 1A. It was used for more accurate segmentation of the pore space that is then later used as a mask for segmenting the fluid phases during the flow experiment. The respective greyscale image from the flow experiment is shown in Figure 1B, where the doped water phase appears in bright and the (undoped) oil in dark. Image processing (filtering and segmentation, see Schlüter, 2010 and Schlüter, 2014) was performed with the Avizo Fire 8.0 software package (Visualization Science Group), following a workflow previously described in (Berg, 2013a; Leu, 2014; Berg, 2015). The reconstructed images were first filtered with a non-local means filter and then registered on the dry scan to apply the mask for the pore space of the rock. Then the fluid phases were segmented using a watershed-based segmentation method followed by a cleanup step using morphological operations. An example of the resulting segmented image is displayed in Figure 1C.

Figure 1: |CT (greyscale) cross-section of dry Gildehauser rock, (b) greyscale image from the flow experiment, showing the (doped) water (bright) and oil (dark) in the pore space, and (C) segmented image from the flow experiment using the segmented pore space from the dry image to mask the pore space.

2.3 Pore-Scale Flow Simulations

The pore-scale imbibition experiment corresponds in principle to an unsteady-state experiment (Mason and Morrow, 2013). However, the accessible length scale of the experiment, i.e. the length of the sample, is smaller than the capillary dispersion zone (Rücker, 2015). Therefore, the Darcy-scale methodology like the analytical JBN inversion technique (Johnson, 1959) or history matching with a reservoir simulator (Berg, 2013a) cannot be used to determine the relative permeability from the experiment. The pore-scale fluid distribution was available from the |CT imaging. Single-phase Navier-Stokes flow simulations were performed on the segmented fluid phases obtained from the |CT imaging during the flow experiment using the software package GeoDict (version 2014; Fraunhofer ITWM, Math2Market GmbH, Kaiserslautern, Germany) directly on the voxel grid from the segmented |CT images. Wetting and non-wetting phase flow were simulated separately for each

time step i of the experiment (with respective wetting phase saturation (Sw ) using a flow rate of 1012 m/s and viscosity of 1 mPa s for the wetting phase and 49.8 mPa s for the non-wetting phase, with a corresponding Reynolds number Re<<1. The resulting permeabilities Kw(Sw ) and Ko(Sw ) were normalized by the absolute permeability k , obtained from a Navier-Stokes flow simulation on the

pore space of the dry sample, yielding (in the connected pathway) relative permeabilities as kr,w(SJ) = Kw(SJ)/Kabs and kr,o(Sj) = Ko(Sj)/Kabs for the wetting and non-wetting phases,

respectively. The obtained relative permeabilities represent the connected pathway flow contribution, but are based on a saturation distribution obtained under dynamic conditions. They do not consider the flux contribution of non-wetting phase ganglion dynamics.

Capillary pressure vs saturation relationships were approximated with a morphological approach (Hilpert and Miller, 2001; Ahrenholz, 2008; Silin, 2010) implemented in GeoDict. Two approximations were applied, which were (i) the morphological approach approximating drainage or imbibition as a sequence of (Young-Laplace) equilibrium states, as opposed to dynamic approaches involving moving liquid-liquid interfaces, and (ii) the Young-Laplace equilibrium states found by a series of morphological erosion and dilation operations using a spherical structuring element. The drainage algorithm starts with a completely water-saturated porous medium and begins with a maximal pore radius R, which is then decreased step by step. In each step, a pore is drained if (i) the pore radius is larger than R, (ii) the entering non-wetting phase is connected by pores larger than R to the nonwetting phase reservoir, and (iii) the replaced wetting phase is connected to the wetting phase reservoir. Conditions (i) and (ii) describe the original algorithm as defined by (Hilpert and Miller, 2001) and used in (Becker, 2008; Vogel, 2005). Condition (iii) introduces a residual wetting phase and was added by (Ahrenholz, 2008). The imbibition algorithm in turn starts with a completely oil-saturated porous medium and minimal pore radius R, which is then increased step by step. In each step a pore is imbibed if (i) the pore radius is smaller than R, (ii) the entering wetting phase is connected to the wetting phase reservoir, and (iii) the replaced non-wetting phase is connected to the non-wetting phase reservoir (Ahrenholz, 2008). This approach is computationally less costly than the (capillary-dominated quasi-static) level-set method (Jettestuen, 2013) and also more efficient than simulating the Navier-Stokes equations for two-phase flow with mobile liquid-liquid interfaces. This approach was first used to obtain a simulated mercury-air (Hg-air) intrusion capillary pressure for comparison with experimentally obtained Hg-air intrusion data (often also referred to as mercury intrusion capillary porosimetry, MICP). In a second step the morphological method was also used to simulate drainage and imbibition. From the resulting pore-scale saturation distribution connected pathway flow relative permeabilities were also obtained in a similar way, as described above. In GeoDict this workflow is automated in the SatuDict module. As evident from approximations taken in the workflow the resulting relative permeabilities were entirely quasi-static.

ig relative

2.4 5teady-5tate Relative Permeability Measurements

For validation of the pore-scale modelling we compare the pore-scale simulation data from the beamline experiment to steady-state relative permeability measurements conducted on a larger dimension sub-sample from the same block. Darcy-scale relative permeability was determined by a steady-state method (Dake, 1978) using a custom-built X-ray saturation measurement setup at Shell in the Netherlands (Kokkedee, 1996; Berg, 2008). Cylindrical rock samples of 5 cm height and 3.8 cm (1.5") diameter were first cleaned using Soxhlett extraction with toluene and, subsequently, an azeotrope mixture of chloroform (84%), methanol (14%), and water (2%). After drying, the porosity $ was determined by helium porosimetry. The absolute permeability to brine Kabs was determined in a flow experiment from the pressure vs flux relationship for multiple flow rates according to eq. 1:

kr,i ■ Kabs APi

q =—---r (1)

Then the sample was mounted inside a sleeve in an X-ray transparent Hassler-type core holder and transferred into the flow setup for steady-state relative permeability measurement by co-injection of saline brine (doped with potassium iodide for increasing the X-ray contrast) and n-decane. The

measurements were conducted at constant flow rate, starting with primary drainage where the fractional flow fw (which is the ratio of water flux over total flux, see also eq. 3) was systematically changed from 100% brine fw = 1) to 100% oil fw = 0), followed by 1st imbibition where the fractional flow fw was systematically changed from 100% oil fw = 0) to 100% brine fw = 1) in 10 saturation steps each. At each saturation step fw was kept constant and the saturation and phase pressures were recorded when steady-state was reached (constant phase pressures and saturations, i.e. in practice variation of the respective averages less than 1% per hour). Pressures were recorded with differential pressure transducers (Kokkedee, 1996). In addition to flow rates and phase pressures, the saturation profiles at 30 positions along the rock sample were also measured, by X-ray transmission measurements (which are calibrated to both 100% brine and oil saturation). The basic procedure and measurement cycles are sketched in Figure 2. For details of the equipment, operating procedures and accuracy we refer to (Kokkedee, 1996; Berg, 2008).

From the resulting phase fluxes qu pressure drops Ap (using pressure data over length L of the central part of the sample only) and saturations, the relative permeability could in principle be directly determined using the two-phase extension of Darcy's law (equation 1; Dake, 1978). However, in order to account for potential capillary end-effects (Huang, 1998), all relative permeability data was obtained by numerical simulation and history matching using Shell's in-house reservoir simulator MoReS (Kokkedee, 1989; Regtien, 1995). More details on the experimental methodology are given in (Kokkedee, 1996; Berg, 2008).

Figure 2: Steady-state relative permeability workflow (A) and schematics of flow setup (B).

2.5 Mercury-Air Intrusion Capillary Pressure Measurements

Mercury porosimetry was used to determine the mercury-air (Hg-air) capillary pressure pc (Stegemeier, 1977). Measurements were performed on smaller sub-samples from the same Gildehauser block using an automatic pore injectivity apparatus (Autopore 9520, Micromeritics). A packing or closure correction is applied during the mercury intrusion, which accounts for packing of mercury around the sample surface prior to entrance of into the pores. The porosity of the sample was determined from the volume of mercury intruded into the sample at 4100 bar, assuming that at this pressure the pore space is completely saturated with mercury.

3 Results and Discussion

3.1 Relative Permeability from ¡uCT Flow Experiment

The basis for the relative permeability calculation from the |CT flow experiment are the (pore-scale) wetting and non-wetting phase saturation distributions obtained at each time step during the imbibition experiment. The respective connectivity of the wetting and non-wetting phases was determined using the Avizo software package. In our experiment both the wetting and non-wetting phases maintain connectivity between inlet and outlet, i.e. both are percolating. However, as the wetting phase saturation increases, the non-wetting phase becomes increasingly disconnected, which is visualized in Figure 3A. Initially the non-wetting phase is mostly connected (red). Then during imbibition the fraction of connected non-wetting phase (red) decreases and the fraction of disconnected non-wetting phase increases (yellow). The decrease in non-wetting phase saturation, loss of internal connectivity and its increasing level of fragmentation into disconnected clusters (Rücker, 2015) causes the permeability of the oil phase to decrease as shown in Figure 3C. As the non-wetting phase saturation decreases, in turn the wetting phase saturation increases and consequently the wetting phase permeability increases, as displayed in Figure 3C. By normalizing the wetting and non-wetting phase connected pathway permeability to the absolute permeability Kabs, the (connected pathway) relative permeabilities were computed and plotted against wetting phase saturation in Figure 4. Relative permeability-saturation functions often show power law-like behavior (Corey, 1954), with values being significant (experimental accuracy on flux measurement is about 103) and relevant for practical field applications over three orders of magnitude. Therefore it is customary in the special core analysis community, which in addition to the hydrology community is also addressed in this paper, to always display relative permeability (and also compare different datasets, validate fits etc.) both on linear and logarithmic scales as is done in Figure 4 and also later figures.

time (s)

Figure 3: (A) Visualization of the non-wetting phase (oil) distribution from |CT images during the flow experiment indicating the contribution of connected (red) and non-connected (yellow) fraction. (B) Saturation of the wetting phase (water) and non-wetting phase. (C) Absolute permeability Kabs, (connected pathway) for oil and water computed from the pore-scale saturation distribution during the dynamic |CT flow experiment.

.Q TO 0 E

-■— k beamline (simulation) # k beamline (simulation) ■■■□■■■■ k steady state imbibition

-O- k steady state imbibition

water saturation S

water saturation S

Figure 4: Comparison between experimental steady-state imbibition relative permeability measurements (open dots)and simulations of data derived from the connected phases of the |CT experiment (closed dots) on both a linear (left) and a logarithmic (right) scale as a function of wetting phase saturation Sw.

Following the main objective for this study, i.e. to investigate to which extent quasi-static methods are sufficient to estimate relative permeability in imbibition, the results from the quasi-static simulations are validated against a reference dataset from steady-state imbibition experiments. In Figure 4 the relative permeabilities from the Darcy-scale steady-state imbibition experiments are superimposed for direct comparison. There is fair agreement between the two datasets in terms of magnitude of relative permeabilities, saturation range, general trends and shape of the curves. The only obvious difference is the saturation endpoint, which for the |CT experiments apparently lies at 10% higher water saturation than for the steady-state experiment.

However, in our experience, such differences in endpoint saturation lie well in the range of sub-sample heterogeneity, and can potentially also be an effect of sample length (Armstrong, 2014). To compensate for differences in endpoint saturation, we replotted in Figure 5 the relative permeability data from Figure 4 in terms of mobile wetting phase saturation

S * =■

S,., - S,..

1 - S_ - S„

representing the saturation interval between the irreducible wetting and non-wetting phase (endpoint) saturations, Sw,c and So,r, respectively. For the |CT experiment, Sw,c = 0.21 and So,r = 0.29, and for the steady-state experiment, Sw, c = 21 and So, r = 0.38 were used, respectively (see also Table 1).

mobile water saturation S* mobile water saturation S*

Figure 5: Comparison of experimental steady-state imbibition relative permeability measurement data with simulations from the connected phases of the |CT experiment on both a linear (left) and a logarithmic (right) scale as a function of the mobile wetting phase saturation S (equation 2).

There is now a better agreement both on the linear and on the logarithmic scale within the experimental error (which for our equipment is approximately 3% for saturation and 10% for relative permeability, see Berg, 2008), except for the water endpoint. It is important to mention here that this comparison is very sensitive to the residual oil saturation So,r in particular for the dataset obtained from the |CT experiment. So,r is very difficult to determine in a reliable way in such dynamic flooding experiments. Centrifuge data would be required for a more reliable determination of So,r, but which was not available.

Apart from the difference in endpoint saturation, i.e. after rescaling to the mobile saturation range, there is good agreement between the steady-state measurement and the quasi-static simulation of the beamline experiment. This means that for this particular case, where the non-wetting phase percolates over the whole saturation range, the contribution of ganglion dynamics (which has been observed in Rucker, 2015) to the overall non-wetting phase flux is minor.

3.2 Impact on Fractional Flow

The question is now whether the difference between the steady-state experimentally measured relative permeability and the quasi-static simulation from the beamline experiment have any significant implications for practical applications such as oil recovery. Without going into the details of fractional flow, a first good indicator for the amount of recoverable oil is the crossover point between krw and kr,o. The crossover point indicates the saturation at which the relative permeability of water and oil are approximately the same and it is a first indicator of the approximate economic limit for oil recovery. In Figure 4 the crossover point for the steady-state data is at Sw = 0.54 whereas for the quasi-static simulation from the beamline experiment the crossover point is at Sw = 0.61, which is significantly different. Using the normalized mobile saturation data in Figure 5, the crossover point for both datasets equals at S = 0.80 (and the exact endpoint saturation So,r can then be determined from a more applicable centrifuge measurement).

For a more detailed analysis, we consider the fractional flow fw for the water phase, which is the ratio of the water flux qw over the total flux (qw + qo). In the limit for vanishing capillary pressure (pc = 0) fw is expressed as (Buckley and Leverett, 1942):

qw + qo 1 + kr o/ Vo

k / li

r ,w r*w

where the water viscosity ¡uw = 1 cP, and the oil viscosity is l = 0.92 cP. In the next step the fractional flow curves fw(Sw) were compared between the quasi-static simulation of the lCT experiment and the experimental steady-state data.

In addition to f directly computed from the individual data points from Figure 4, also smoother f curves were computed from fits of the relative permeability data with Corey functions (Corey, 1954; Masalmeh, 2006):

kr,w (Sw ) = kr,wL "(S* Y

kr,o (Sw ) = kr.

'I1-S„

• (1 - S * f

with the mobile saturation S* from eq. (2). The results of the fit are displayed in Figure 6, indicating a valid representation of the steady-state data by the Corey function. The respective parameters of the Corey fits, i.e. endpoint saturations and endpoint relative permeability k I and k I and

the Corey exponents nw and no for water and oil, respectively, are listed in Table 1.

The Corey function parameterizes relative permeability with a smooth curve that follows a power law functionality and the last data point, i.e. the endpoint, should be consistent with the overall trend of the curve representing the other data points. Assuming that relative permeability is well described by Corey functions, the comparison presented in Figure 6 provides an additional level of consistency and reduces the uncertainty with respect to the endpoint saturation So,r. For the simulation data of the ^CT experiment, the oil relative permeabilities are represented well by Corey functions but the water relative permeabilities at low saturations clearly deviate from a Corey representation (deviation larger than the experimental error). In particular in the saturation range S < 0.6 the agreement for kro is good but for S* > 0.6 the agreement is modest at best, in particular for the water endpoint. The origin of this particular discrepancy is not clear. However, the impact on the fractional flow curve in Figure 7 and the breakthrough saturations and recoveries are within the experimental uncertainty.

Table 1: Corey fit parameters and results of a Buckley-Leverett calculation for the relative data from Figure 4. H

ere Swc is

the connate water saturation, Sor the residual oil saturation, kr,o(Sw,c) the oil relative permeability endpoint (at connate water saturation), kr,w(1-So,r) the water relative permeability endpoint (i.e. at residual oil saturation), and nw and no are the water and oil Corey exponents, respectively. Sw@BT is the water saturation at breakthrough, fw@BT the respective fractional flow, and Recovery@BT the respective oil recovery.

tJw, c So, r kr, o(Sw,c) kr,w(1-So,r) nw no Sw@BT fw@BT Recovery@BT

Simulation of ^CT experiment 0.21 0.29 0.56 0.16 5.5 1.7 0.69 0.97 0.62

Steady-state 0.21 0.38 0.52 0.06 2.4 1.6 0.62 0.99 0.52

"§0.6 E <¡3

<u 0.2

E <¡3

0.0 0.0

water saturation S

Figure 6: Corey fits of the simulated |CT-based relative permeability (top row) and the steady-state imbibition relative permeability (bottom row) on both a linear (left) and a logarithmic (right) permeability axis.

ra c o

1.0 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0.0

beamline

- ... Corey fit BL tangent steady-state

- ... Corey fit BL tangent

0.2 0.4 0.6 0.8

water saturation S

1.0 0.0

0.2 0.4 0.6 0.8

mobile water saturation S*

Figure 7: Comparison of the fractional flow curves for the experimental steady-state imbibition relative permeability (measurement) and simulation from the connected phases of the pore scale as a function of the water saturation Sw (left) and the mobile wetting phase saturation S (right). Here also the position of the shock front is indicated as the point where the tangent meets the fractional flow curve following the Buckley-Leverett (BL) formalism.

In Figure 7, the fractional flow curves of the quasi-static simulation of the ^CT experiment were compared with the experimental steady-state data using both the absolute saturation Sw and the mobile saturation S scales. For 1D displacement in the absence of any capillary pressure, key parameters such as breakthrough water saturation (Sw@BT), water cut at breakthrough (fw@BT) and recovery at breakthrough (Recovery@BT) of the shock front can be computed analytically from the fractional flow curve fw(Sw) using the Buckley-Leverett formalism (Buckley and Leverett, 1942; Dake, 1978). The graphical representation of these calculations is a tangent from Sw,c to the fractional flow curve included in Figure 7. The saturation at which the tangent intersects with fw is the breakthrough saturation, and the respective value of fw is the water cut at breakthrough. The results of the Buckley-Leverett calculation for the steady-state experiment and simulation results on basis of the ^CT experiment are listed in Table 1. Without normalizing to mobile saturation, there is a difference of 10% in saturation units for recovery at breakthrough between the steady-state experiment and the simulation results, which is roughly the difference in endpoint saturation. On the mobile saturation scale S* the fractional flow curves and the respective tangents are very similar. Consequently, in the hypothetical case that the endpoints were identical, the respective breakthrough saturations and recovery at breakthrough would be identical within less than 5% for both the steady-state experiment and simulation results.

3.3 Capillary Pressure from the Dry Pore Space Using a Morphological Approach

A comparison between the experimentally measured (MICP experiment) and simulated (using the morphological approach in GeoDict, Becker, 2008; Hilpert and Miller, 2001; Vogel, 2005) Hg-air capillary pressurepc data is displayed in Figure 8. Note that both cases represent capillary-dominated techniques where the contribution of viscous forces is absent or at least negligible. There is a good agreement between the simulated and the MICP data in a wetting phase saturation range above 0.2, which is similar to observations presented by Hussain, 2014. The discrepancies at lower wetting phase saturations are related to limitations by the imaging resolution, where the pore space is not fully resolved (Leu, 2014).

wetting phase saturation

Figure 8: Simulated Hg-air intrusion capillary pressure with morphological approach using GeoDict and comparison with experimental MICP data for Gildehauser sandstone rock. In addition, the GeoDict oil-water "imbibition" capillary pressure simulation scaled to Hg-air by the interfacial tension ratio between oil-water and Hg-air is plotted for comparison.

3.4 Relative Permeability from Morphological Approach

Apart from end-point saturations, relative permeability obtained from quasi-static simulations based on experimentally measured pore-scale fluid distributions show good agreement with steady-state relative permeability reference measurements. One may argue that this is good news because for practical applications the computationally much cheaper quasi-static simulation could be used for predicting relative permeability. However, this approach is not practical because the pore-scale fluid distributions of wetting and non-wetting phases are required, which for the presented data requires synchrotron-based lCT flow experiments. Therefore, as the next step a morphological approach (Becker, 2008; Hilpert and Miller, 2001; Vogel, 2005) was used to simulate the capillary-dominated invasion process. As already shown in Figure 8, this approach gives excellent agreement between simulated and measured mercury-air intrusion capillary pressure-saturation relationships. Similar observations have been made for Berea sandstone rock in a previous study (Leu, 2014). However, mercury-air intrusion is a drainage process. The important question is now how well this approach works for imbibition. Therefore in the next step we simulate pore-scale fluid distributions using the (capillary-controlled) morphological approach (Becker, 2008; Hilpert and Miller, 2001; Vogel, 2005) in "imbibition mode" as implemented in GeoDict. Then we compute the quasi-static hydraulic conductivities on the resulting pore-scale fluid distribution analogous to the approach taken with the LCT data. The final results are compared to the steady-state relative permeability reference data.

The results show that the morphological approach using the "imbibition mode" setting provide relative permeability curves that are more representative of a drainage process than an imbibition process. As displayed in Figure 9, the resulting imbibition relative permeability shows very large discrepancies from the measured relative permeability and the simulation from the beamline data.

Figure 9: Comparison of experimental steady-state imbibition relative permeability (measurement) and simulation from the connected phases of the pore-scale experiment to the morphological approach on a linear (left) and a logarithmic (right) scale as a function of the mobile wetting phase saturation S .

However, the water relative permeability from the morphological approach matches well with the drainage steady-state relative permeability as shown in Figure 10. There are still discrepancies with the oil relative permeability, which indicates that the morphological approach has fundamental difficulties in "simulating" connected oil phases, in particular at low oil saturations.

0.9 0.8 0.7

£ 0.6

§ 0.5 <D

<D 0.4 >

13 0.3

■■■□■■■■ k steady state drainage -O- k steady state drainage

k morphological approach -A- k morphological approach

0.01 -

mobile water saturation S*

obile water saturation S*

Figure 10: Comparison of experimental steady-state drainage relative permeability (measurement) and the morphological approach on a linear (left) and logarithmic scale (right) as a function of the mobile wetting phase saturation S .

Also, the capillary pressure data from Figure 8, where the scaled oil-water imbibition capillary pressure superimposes with the Hg-air drainage capillary pressure, indicates that the morphological approach approximates more a drainage process with respect to the computed pore-scale fluid distribution. The reason why capillary pressure is predicted reasonably well by the quasi-static method is that the morphological method describes correctly in which pore the non-wetting phase is distributed in individual pores, also probing directly connected pore throats. Also, in a drainage process the conductivity is captured sufficiently well because it is an invasion process where, by design of the algorithm used, both phases are connected.

For imbibition, however, the long-range connectivity of the non-wetting phase is not described very well because processes like snap-off and corner flow are not captured. In reality, there is mass exchange between connected and disconnected fluid phases (Hilfer, 1998; Hilfer, 2006a; Hilfer, 2006b), which is not only adding to the net flux but also influences the time evolution of the pore-scale fluid distributions, including the configuration of the connected phases. The morphological approach used here employs local physics to determine the pore-scale saturation distribution (and hence is conceptually similar to quasi-static pore network modelling which uses local rules). Approaches based on local physics cannot capture the non-local and cooperative processes that have been recently observed in drainage (Armstrong, 2013; Armstrong, 2015; Andrew, 2015) and imbibition (Rücker, 2015; Holtzman, 2015). However, our findings show that the local approximation has a more dramatic effect on imbibition than on drainage. Our work clearly shows that this has less to do with the quasi-static computation of the flow field based on a given saturation distribution but rather the approximation of the saturation distribution by a quasi-static morphological approach. A possible explanation can be offered by considering the different pore-scale displacement processes for drainage and imbibition. While for drainage frontal displacement and to some extent cooperative filling are the most relevant process, for imbibition the snap-off and corner flow processes are fundamental to what ultimately determines the connectivity of the non-wetting phase. Snap-off is a

process that is influenced by both the local pore and fluid geometry in pore throats, but also the nonlocal, long-range connectivity of the non-wetting phase which is required to provide the volume of wetting fluid to allow the snap-off to occur (Mohanty, 1987). By applying only local rules, this longrange connectivity is not considered and only based on local rules snap-off processes are permitted even if the long-range connectivity may suppress them. Also, to capture the physics of corner the flow the tiniest crevices of the pore-structure must be imaged and simulated with sufficient resolution. Different realizations of the pore space corners would result in wetting film advancement along different paths and thus ultimately influence the sequence of pore filling events. Therefore the resulting pore-scale fluid distribution is different (and difficult to capture for imbibition), which ultimately leads to a significantly different connectivity and hence relative permeability of particularly the non-wetting phase.

With respect to the data from the morphological approach it is also important to note that the oil relative permeability is only accessible until saturation values of 5* « 0.4. This is because for larger water saturations, i.e. lower oil saturations, the oil phase is not percolating any more, as also illustrated in Figure 11, and thus the connected pathway hydraulic conductivity is zero. In the morphological approach the saturation value of 5w « 0.36-38 or 5* « 0.4 corresponds to a percolation threshold for the non-wetting phase whereas in the experiment the non-wetting phase percolates up to much higher saturations of 5w « 0.70.

nw phase percolates

connected nw-phase disconnected nw-phase

nw phase percolates

Figure 11: Comparison of the pore-scale fluid distributions from the beamline experiment (top) with the morphological approach from GeoDict in imbibition (bottom) as a function of the wetting phase saturation Sw. Red indicates a connected and yellow a disconnected non-wetting phase. For the morphological approach the non-wetting phase disconnects i.e. stops percolating at much lower water saturations (around 36-38%) compared to the beamline data where the non-wetting phase percolates up to 70% in this particular experiment.

Currently, it is not fully clear whether it can be expected in all cases that connected pathway flow will dominate over ganglion dynamics and that the percolation threshold for the non-wetting phase ranges up to 5w « 0.70. Recent studies by Reynolds et al. (Reynolds, 2015) provide an example with lower percolation threshold and hence more contribution of ganglion dynamics. There the connected pathway flow approach would not succeed to compute relative permeability beyond the percolation threshold.

4 Summary and Conclusions

The ability to predict parameters for Darcy-scale flow from quasi-static pore-scale flow simulations conducted directly on the pore space imaged by ^CT was investigated. Previous studies have shown reasonably good agreement with Darcy-scale measurements for both drainage and imbibition when using a capillary-viscous flow simulation, which is computationally very expensive (Vogel, 2005; Ramstad, 2009; Ferrari, 2013; Koroteev, 2014). Therefore, the focus here is on quasi-static simulation, which comes at a substantially lower computational cost.

Previous studies (Hussain, 2014) have demonstrated that for drainage good agreement with Darcy-scale measurements can be obtained using a morphological approach to estimate the pore-scale fluid distribution and then conducting quasi-static flow simulations on the connected fluid pathways. We find that for imbibition this approach works well if the pore-scale fluid distributions are available from a flow experiment with ^CT pore-scale imaging. Apart from a 10% difference in the endpoint saturation between twin samples, which can also have different origins than the applied methodology, the relative permeability computed on the connected pathway from a synchrotron beamline flow experiment matches closely with the Darcy-scale steady-state measurement. This means that for a given pore-scale fluid distribution where the non-wetting phase percolates over the whole saturation range, the contribution of ganglion dynamics (Rucker, 2015) to the overall flux is very minor. However, it is not fully clear whether it can be expected in all cases that connected pathway flow will dominate over ganglion dynamics. Recent studies by Reynolds et al. (Reynolds, 2015) provide an example that seems to be more dominated by ganglion dynamics and there at least the connected pathway flow approach would not succeed to compute relative permeability, at least not at lower oil saturations.

In summary, we find that quasi-static flow simulations give a good approximation of relative permeability when the pore-scale fluid distributions at each saturation step are known and the oil is mainly connected. If the pore-scale fluid distributions are not available, the approximation with a morphological approach does provide a good basis for the computation of drainage capillary pressure (which matches with measured Hg-air intrusion data very well) but not for imbibition relative permeability. The reason is that quasi-static morphological approaches are not capable of predicting such pore-scale saturation distributions for imbibition. The relative permeability computed from the respective pore-scale fluid distribution using the connected pathway conductivity of the two fluid phases is much closer to a measured drainage relative permeability.

Ton Blok is gratefully acknowledged for the MICP measurements. Sebastiaan Pieterse and Niels Brussee are gratefully acknowledged for the steady-state relative permeability measurements. Fons Marcelis and Ab Coorn are acknowledged for sample preparation and Soxhlett cleaning. We acknowledge the Paul Scherrer Institut, Villigen, Switzerland for provision of synchrotron radiation beamtime at the TOMCAT beamline of the SLS and would like to thank Sally Irvine and Marco Stampanoni for assistance.

6 References

Ahrenholz, B., Tölke, J, Lehmann, P., Peters, A., Kaestner, A., Krafczyk, M. and Durner, W., Prediction of capillary hysteresis in a porous material using lattice-Boltzmann methods and comparison to experimental data and a morphological pore network model. Adv. Water Resour. 31, 1151-1173 (2008)

Al Ibrahim, M. A., Hurley, N. F., Zhao, W., Acero-Allard, D., An Automated Petrographic Image Analysis System: Capillary Pressure Curves Using Confocal Microscopy, SPE Annual Technical Conference and Exhibition held in San Antonio, Texas, USA, 8-10 October 2012, SPE 159180 (2012).

Andrew, M., Menke, H., Blunt, M. J., Bijeljic, B., The Imaging of Dynamic Multiphase Fluid Flow Using Synchrotron-Based X-ray Microtomography at Reservoir Conditions, Transport in Porous Media 110, 1-24 (2015).

Armstrong, R. T. and Berg, S., Interfacial velocities and capillary pressure gradients during Haines jumps. Phys. Rev. E 88, (2013)

Armstrong, R. T., Georgiadis, H., Ott, H. Klemin, D., Berg, S., Pore-scale mobilization of non-wetting phase ganglia and critical capillary number for onset of macroscopic capillary desaturation, Geophys. Res. Lett. 41, 16, 2014.

Armstrong, R. T., Evseev, N., Koroteev, D., Berg, S., Velocity field during Haines jumps in porous media, Advances in Water Resources, DOI:10.1016/j.advwatres.2015.01.008 (2015).

Avraam, D. G., Payatakes, A. C., Flow regimes and relative permeabilities during steady-state two-phase flow in porous media, Journal of Fluid Mechanics 293, 207-236 (1995).

Becker, J. Schulz, V. Wiegmann, A., Numerical determination of two-phase material parameters of a gas diffusion layer using tomography images. J. Fuel Cell Sci. Tech. 5, (2008)

Berg, S., Cense, A. W., Hofman, J. P., Smits, R. M. M., Two-Phase Flow in Porous Media with Slip Boundary Condition, Transport in Porous Media 74(3), 275-292 (2008).

Berg, S. Ott, H. Klapp, S. A. Schwing, A. Neiteler, R. Brussee, N. Makurat, A. Leu, L. Enzmann, F. Schwarz, JO. Kersten, M. Irvine, S. Stampanoni, M., Real-time 3D imaging of Haines jumps in porous media flow. In PNAS 10, 3755-3759 (2013a)

Berg, S., Oedai, S., Ott, H., Displacement and Mass Transfer Between Saturated and Unsaturated CO2-Brine Systems in Sandstone, International Journal of Greenhouse Gas Control 12, 478-492, (2013b).

Berg, S. Armstrong, R. Ott, H. Georgiadis, A. H. Klapp, S. A. Schwing, A. Neiteler, R. Brussee, N. Makurat, A. Leu, L. Enzmann, F. Schwarz, J-O. Wolf, M. Khan, F. Kersten, M. Irvine, S. Stampanoni, M., Multiphase flow in porous rock imaged under dynamic flow conditions with fast x-ray computed microtomography. Petrophysics 55(4), 304-312 (2014).

Berg, S., Rücker, M., Amstrong, R. T., Georgiadis, A., Ott, H., Leu, L., Wolf, M., Khan, F., Enzmann, F., Kersten, M., Onset of Oil Mobilization and Non-wetting Phase Cluster Size Distribution, Petrophysics 56(1), 15-22 (2015).

Bijeljic, B., Raeini, A., Mostaghimi, P., Blunt, M. J., Predictions of non-Fickian solute transport in different classes of porous media using direct simulation on pore-scale images. Phys. Rev. E 87, (2013).

Blunt, M. J., Jackson, M. D., Piri, M., Valvatne, P. H., Detailed physics, predictive capabilities and macroscopic consequences for pore-network models of multiphase flow, Advances in Water Resources 25, 1069-1089 (2002).

Blunt M. J. Bijeljic, B. Dong, H. Gharbi O. Iglauer, S. Mostaghimi, P. Paluszny, A. Pentland, C., Pore-scale imaging and modeling. Adv. Water Resour. 51, 197-216 (2013).

Boek, E. S., Zacharoudiou, I., Gray, F., Shah, S. M., Crawshaw, J. P., Yang, J., Multiphase Flow and Reactive Transport at the Pore Scale Using Lattice-Boltzmann Computer Simulations, SPE Annual Technial Conference and Exhibition held in Amsterdam, The Netherlands, 27-29 October 2014, SPE-170941, 2014.

Buckley, S. E., Leverett, M. C., Mechanism of Fluid Displacement in Sands, Transactions of the AIME 146(1), 107-116 (1942).

Chukwodozie, C., Tyagi, M., Pore Scale Inertial Flow Simulations in 3-D Smooth and Rough Sphere Packs using Lattice Boltzmann Method, AIChE Journal 59(12), 4858-4870 (2013).

Coles, M. E., Hazlett, R. D., Spanne, P., Soll, W. E., Muegge, E. L., Jones, K. W., Pore level imaging of fluid transport using synchrotron X-ray microtomography. In J. Pet. Sci. Eng. 19, 55-63 (1998).

Corey, A.T., The Interrelation Between Gas and Oil Relative Permeabilities. Prod. Monthly 19 (1), 38-41 (1954).

Dake, L. P., Fundamentals of Reservoir Engineering. Elsevier B.V. (1978).

Datta, S. S., Dupin, J.-B., Weitz, D. A., Fluid breakup during simultaneous two-phase flow through a three-dimensional porous medium, Physics of Fluids 26, 062004 (2014).

Demianov, A., Dinariev, O., Evseev, N., Density Functional Modelling in Multiphase Compositional Hydrodynamics, The Canadian Journal of Chemical Engineering 89, 206-226 (2011).

Dubelaar, C. W., Nijland, T. G., Early Cretaceous Obernkirchen and Bentheim Sandstones from Germany used as dimension stone in the Netherlands: geology, physical properties, architectural use and comparative weathering. Geological Society, London, Special Publications, 416, SP416-13, (2015).

Dullien, F. A.L., Porous Media: Fluid Transport and Pore Structure, Academic Press, New York City (1979).

Ferrari, A., Lunati, I., Direct numerical simulations of interface dynamics to link capillary pressure and total surface energy, Advances in Water Resources 57, 19-31 (2013).

Fusseis, F., Steeb, H., Xiao, X., Zhu, W.-I., Butler, I. B., Elphick, S., Mader, U., A low-cost X-ray-transparent experimental cell for synchrotron-based X-ray microtomography studies under geological reservoir conditions,

Journal of Synchrotron Radiation 21, 251-253 (2014).

Georgiadis, A., Berg, S., Maitland, G. Makurat, A., Ott, H., Pore-scale micro-CT Imaging: non-wetting phase cluster size distribution during drainage and imbibition. Phys. Rev. E 88, (2013).

Ghanbarzadeh, S., Prodanovic, M., Hesse, M., Percolation and Grain Boundary Wetting in Anisotropic Textually Equilibrated Pore Networks, Physical Review Letters 113, 048001 (2014).

Gray, W. G., Dye, A. L., McClure, J. E., Pyrak-Nolte, L. J., Miller, C. T., On the dynamics of two-fluid flow in porous media, Water Resources Research, DOI 10.1002/2015WR016921, 2015.

Hammond, P. S., Unsal, E., A Dynamic Pore Network Model for Oil Displacement by Wettability-Altering Surfactant Solution, Transport in Porous Media 92(3), 789-817 (2012).

Hilfer, R., Macroscopic equations of motion for two-phase flow in porous media, Phys. Rev. E: Stat. Nonlinear Soft Matter Phys. 58(2), 2090-2096 (1998).

Hilfer, R., Macroscopic capillarity and hysteresis for flow in porous media, Phys. Rev. E: Stat. Nonlinear Soft Matter Phys. 73, 016307 (2006a).

Hilfer, R., Macroscopic capillarity without a constitutive capillary pressure function, Physica A 371, 209-225 (2006b).

Hilfer, R., and P. E. 0ren, Dimensional analysis of pore scale and field scale immiscible displacement, Transp. Porous Media 22(1), 53-72 (1996).

Hilpert, M. and Miller, T. C., Pore-morphology-based simulation of drainage in totally wetting porous media. In

Adv. Water Resour. 24, 243-255 (2001).

Holtzman, R., Segre, E., Wettability Stabilizes Fluid Invasion into Porous Media via Nonlocal, Cooperative Pore Filling, Physical Review Letters 115, 164501 (2015).

Huang, D. D., Honarpour, M. H., Capillary end effects in coreflood calculations, Journal of Petroleum Science and Engineering 19, 103-117 (1998).

Hussain, F., Pinczweski, W. V., Cinar, Y., Arns, J. Y., Arns, C. H., Turner, M. L., Computation of Relative Permeability from Imaged Fluid Distributions at the Pore Scale, Transp. Porous Med. 104, 91-107 (2014).

Iglauer, S. Paluszny, a. Pentland, C. Blunt, M. J., Residual CO2 imaged with X-ray micro-tomography. Geophys. Res. Lett. 38, (2011).

Iglauer, S., Ferno, M. A., Shearing, P., Blunt, M. J., Comparison of residual oil cluster size distribution, morphology and saturation in oil-wet and water-wet sandstone, Journal of Colloid and Interface Science 375, 187-192 (2012).

Jacqmin, D., Calculation of Two-Phase Navier-Stokes Flows Using Phase-Field Modeling, Journal of Computational Physics 155, 96-127 (1999).

Jettestuen, E., Helland, J. O. and Prodanovic, M., A level set method for simulating capillary-controlled displacements at the pore scale with nonzero contact angles, Water Resources Research 49, 1-17 (2013).

Joekar-Niasar, V., Hassanizadeh, S. M., Uniqueness of Specific Interfacial Area-Capillary Pressure-Saturation Relationship Under Non-Equilibrium Conditions in Two-Phase Porous Media Flow, Transport in Porous Media 94, 465-486 (2012).

Johnson, E. F., Bossler, D.P., and Naumann, V.O., Calculation of Relative Permeability from Displacement Experiments, Trans. AIME, 216, 370-372, 1959.

Kokkedee, J.A., Boom, W., Frens, A.M., Maas, J.G., Improved special core analysis: scope for a reduced residual oil saturation. Society of Core Analysis Conference Paper. 9601, 1-13 (1996).

Koroteev, D., Dinariev, O., Evseev, N., Klemin, D., Nadeev, A., Safonov, S., Gurpinar, O., Berg, S., van Kruijsdijk, C., Armstrong, R. T., Myers, M. T., Hathon, L., de Jong H., Direct Hydrodynamic Simulation of Multiphase Flow in Porous Rock, Petrophysics 55(4), 294-303, 2014.

Lake, L. W., Enhanced Oil Recovery, Prentice Hall, 1989.

Landry, C. J., Karpyn, Z. T., Piri, M., Pore-scale analysis of trapped immiscible fluid structures and flui interfacial areas in oil-wet and water-wet bead packs, Geofluids 11, 209-227 (2011).

Landry, C. J., Karpyn, Z. T., Ayala, O., Pore-Scale Lattice Boltzmann Modeling and 4D X-ray Computed Microtomography Imaging of Fracture-Matrix Fluid Transfer, Transport in Porous Media 103, 449-468 (2014).

Lenormand, R., Touboul, E., Zarcone, C., Numerical models and experiments on immiscible displacements in porous media, Journal of FluidMech. 189, 165-187 (1988).

Leu, L., Berg, S., Enzmann, F., Armstrong, R. T., Kersten, M., Fast X-ray micro-tomography of multi-phase flow in Berea sandstone: a sensitivity study on image processing, Transport in Porous Media 105 (2) 451-469, (2014).

Leverett, C. M., Capillary behaviour in porous solids. In Transact. AIME 142, (1941)

Manuel, D., Van Loo, D., Masschaele, B., Van den Bulcke, J., Van Acker, J., Cnudde, V., Van Hoorebeke, L., Recent micro-CT scanner developments at UGCT, Nuclear Instruments and Methods in Physics Research B 324, 35-40 (2014).

Masalmeh, S.K. and Jing, X.D., Capillary Pressure Characteristics of Carbonate, Reservoirs: Relationship between Drainage and Imbibition Curves, SCA 2006-16, Reviewed Proc. International Symposium of the Society of Core Analysts, Trondheim, Norway, 12-16 September, (2006).

Mason, G., Mobilization of Oil Blobs in the Pore Space of Random Sphere Packings, Chemical Engineering Science 38(9), 1455-1460 (1983).

Mason, J. and Morrow, N. R., Developments in spontaneous imbibition and possibilities for future work, Journal of Petroleum Science and Engineering 110, 268-293 (2013).

Mohanty, K. K., Davis, H. T., Scriven L. E., Physics of oil entrapment in water-wet rock. SPE Res Eval & Eng 2(1):113-128 (1987).

Mostaghimi, P. Blunt, M. J. Bijeljic, B., Computations of absolute permeability on micro-CT images. Mathem. Geosci. 45, 103-125 (2013).

Niessner, J., Berg, S., Hassanizadeh, S. M., Comparison of two-phase Darcy's law with a thermodynamically consistent approach, Transport in Porous Media 88(1), 133-148 (2011)

Porter, M. L., Schaap, M. G., Wildenschild, D., Lattice-Boltzmann simulations of the capillary pressure-saturation-interfacial area relationship for porous media, Advances in Water Resources 32(11), 1632-1640, 2009.

Prodanovic, M., Lindquist, W. B., Serigh, R. S., Porous structure and fluid partitioning in polyethylene cores from 3D X -ray microtomographic imaging, Journal of Colloid and Interface Science 298, 282-297 (2006).

Raeesi, B., Morrow, N., Mason, G., Pore Network Modeling of Experimental Capillary Pressure Hysteresis Relationships, International Symposium of the Society of Core Analysts held in Napa Valley, California,, USA, 16-19 September, (2013).

Raeini, A. Q., Bijeljic, B., Blunt, M. J., Numerical Modelling of Sub-pore Scale Events in Two-Phase Flow Through Porous Rock, Trans. Porous Med. 101, 191-213 (2014a).

Raeini, A. Q., Blunt, M. J., Bijeljic, B., Direct simulations of two-phase flow on micro-CT images of porous media and upscaling of pore-scale forces, Advances in Water Resources 74, 116-126 (2014b).

Ramstad, T., Oren, P.-E., Bakke, S., Simulation of Two Phase Flow in Reservoir Rocks Using a Lattice Boltzmann Method, SPE 124617, SPE Annual Technical Conference and Exhibition held in New Orleans, Lousiana, USA, 4-7 October 2009.

Raoof, A., Nick, H. M., Hassanizadeh, S. M., Spiers, C. J., PoreFlow: A complex pore-network model for simulation of reactive transport in variably saturated porous media, Computers & Geosciences 61, 160-174 (2013).

Regtien, J. M. M., Por, G. J. A., van Stiphout, M. T. and van der Vlugt, F. F., Interactive reservoir simulation,

Proceedings of the 13th SPE Symposium on Reservoir Simulation, San Antonio, TX, 1995, pp. 545-552, SPE-

Reynolds, C., Menke, H., Andrew, M., Blunt, M., Krevor, S., Observations of connected and disconnected nonwetting phase flow during co-injection of N2 and brine, Proceedings of the 7th International Conference on Porous Media & Annual Meeting held in Padova, Italy, May 18-21, 2015.

Rücker, M., Berg, S., Armstrong, R. T., Georgiadis, A., Ott, H. Schwing, A., Neiteler, R., Brussee, N., Makurat, A., Leu, L. Wolf, M., Khan, F., Enzmann, F., Kersten, M., From Connected Pathway Flow to Ganglion Dynamics, Geophysical Research Letters 42, doi:10.1002/2015GL064007 (2015).

Schlüter, S. Weller, U. Vogel, H.-J., Segmentation of X-ray microtomography images of soil using gradient masks. Comp. Geosci 36, 1246-51 (2010)

Schlüter, S. Sheppard, A. Brown, K. Wildenschild, D., Image processing of multiphase images obtained via X-ray microtomography: A review. Water Resour. Res. 50, 3615-3639 (2014)

Silin, D., Tomutsa, L. Benson, S. Patzek, T. W., Microtomography and pore-scale modeling of two-phase fluid distribution. Transp. Porous Media 86, 495-515 (2010).

Stegemeier, G. L., Mechanisms of Entrapment and Mobilization of Oil in Porous Media, in Improved Oil Recovery by Surfactant and Polymer Flooding, Shah, D. O. and Schechter, R. S. ed., Academic Press, New York (1977).

Turner, M., Knüfing, L., Arns, C., Sakellariou, A., Senden, T., Sheppard, A., Sok, R., Limaye, A., Pinczewski, W.V., Knackstedt, M., Three-dimensional imaging of multiphase flow in porous media. Phys. A 339, 166-172 (2004)

Vogel, H. J. Tölke, J. Schulz, V. P. Krafczyk, M. Roth, K., Comparison of a lattice-Boltzmann model, a full morphology model, and a pore network model for determining capillary pressure-saturation relationships.

Vadose Zone J. 4 , 380-388 (2005).

Wildenschild, D., Sheppard, A. P., X-ray imaging and analysis techniques for quantifying pore-scale structure and processes in subsurface porous medium systems, Advances in Water Resources 51, 217-246 (2013).

Yadav, G. D., Mason, G., The Onset of Blob Motion in a Random Sphere Packing Caused by Flow of the Surrounding Fluid, Chemical Engineering Science 38(9), 1461-1465 (1983).

29146.

Deschamps, H., Oughanem, R., Maire, E., Mokso, R., Oil ganglia dynamics in surfactant flooding captured by ultra-fast x-ray microtomography, International Core Analysts held in Avignon, France, 11-18 September, 2014, SCA2014-023

A Representative Elementary Volume

A representative elementary volume (REV) analysis with respect to porosity and permeability for the Gildehauser sandstone sample used for the synchrotron beamline flow experiments was performed using the ^CT image from Figure 1 following the methodology by (Al Ibrahim, 2012). For increasing window size, porosity and permeability values were determined and plotted as a function of window size as shown in Figure A1. In addition the mean + 1 standard deviation (standard deviation is obtained from multiple sub-samples at same window size) is plotted. According to (Al Ibrahim, 2012) the intersection between mean and mean + 1 standard deviation is an estimate for the REV. The results in Figure A1 show that the sample size of 4 mm diameter is larger than the (single phase) REV. Note that while for porosity the mean becomes independent of the window size for windows sizes larger than the REV, for permeability this is not the case because it is not a linear property. Also note that the 2-phase REV can be substantially larger than the single-phase REV (Georgiadis, 2013

.Q CO CD

.: Repre

! mean F mean + 1 stdev L ---□ — -■

—-——................

-■-■-■—.....'

window size l (^m)

Figure A1: Representative Elementary Volume (REV) analysis with respect to porosity (top panel) and permeability (bottom panel) for Gildehauser sandstone rock following the methodology by (Al Ibrahim, 2012). The REV is estimated from the point where the mean and the mean + 1 standard deviation intersection.

Graphical abstract

-m- krw beamline (simulation)

ko beamline (simulation)

k steady state imbibition r,w

■ o- k o steady state imbibition

0.2 0.4 0.6 0.8 water saturation S