Regular Article
Analytical modeling and correction of steady state relative permeability experiments with capillary end effects – An improved intercept method, scaling and general capillary numbers
^{1}
Department of Energy Resources, University of Stavanger, 4021 Stavanger, Norway
^{2}
The National IOR Centre of Norway, University of Stavanger, 4021 Stavanger, Norway
^{*} Corresponding author: pal.andersen@uis.no
Received:
14
May
2021
Accepted:
29
July
2021
Steady state relative permeability experiments are performed by coinjection of two fluids through core plug samples. Effective relative permeabilities can be calculated from the stabilized pressure drop using Darcy’s law and linked to the corresponding average saturation of the core. These estimated relative permeability points will be accurate only if capillary end effects and transient effects are negligible. This work presents general analytical solutions for calculation of spatial saturation and pressure gradient profiles, average saturation, pressure drop and relative permeabilities for a core at steady state when capillary end effects are significant. We derive an intuitive and general “intercept” method for correcting steady state relative permeability measurements for capillary end effects: plotting average saturation and inverse effective relative permeability (of each phase) against inverse total rate will give linear trends at high total rates and result in corrected relative permeability points when extrapolated to zero inverse total rate (infinite rate). We derive a formal proof and generalization of the method proposed by Gupta and Maloney (2016) [SPE Reserv. Eval. Eng. 19, 02, 316–330], also extending the information obtained from the analysis, especially allowing to calculate capillary pressure. It is shown how the slopes of the lines are related to the saturation functions allowing to scale all test data for all conditions to the same straight lines. Two dimensionless numbers are obtained that directly express how much the average saturation is changed and the effective relative permeabilities are reduced compared to values unaffected by end effects. The numbers thus quantitatively and intuitively express the influence of end effects. A third dimensionless number is derived providing a universal criterion for when the intercept method is valid, directly stating that the end effect profile has reached the inlet. All the dimensionless numbers contain a part depending only on saturation functions, injected flow fraction and viscosity ratio and a second part containing constant known fluid, rock and system parameters such as core length, porosity, interfacial tension, total rate, etc. The former parameters determine the saturation range and shape of the saturation profile, while the latter number determines how much the profile is compressed towards the outlet. End effects cause the saturation profile and average saturation to shift towards the saturation where capillary pressure is zero and the effective relative permeabilities to be reduced compared to the true relative permeabilities. This shift is greater at low total rate and gives a false impression of ratedependent relative permeabilities. The method is demonstrated with multiple examples. Methodologies for deriving relative permeability and capillary pressure systematically and consistently, even based on combining data from tests with different fluid and core properties, are presented and demonstrated on two datasets from the literature. The findings of this work are relevant to accurately estimate relative permeabilities in steady state experiments, relative permeability end points and critical saturations during flooding or the impact of injection chemicals on mobilizing residual phase.
© Pål Ø. Andersen, published by IFP Energies nouvelles, 2021
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Nomenclature
RomanC_{s}, C_{i}: Saturation and relative permeability slopes in the intercept method, −
F_{i} : Injected phase fraction, −
f_{w} : Water fractional flow function, −
J : Scaled capillary pressure, −
J_{1}, J_{2}, J_{3}: Jfunction coefficients, −
k_{ri} : Phase relative permeability, −
K : Absolute permeability, m^{2}
k_{1}, k_{2}: Jfunction parameters, −
n : Exponent for approximating end effect saturation profile, −
n_{1}, n_{2}: Jfunction exponents, −
n_{i} : Phase Corey exponent, −
n_{i1}, n_{i2}: Phase Corey exponent end values, −
N_{0} : capillary number (viscous to capillary forces), −
P_{c} : Capillary pressure, Pa
: Water saturation at which capillary pressure is zero, −
: Water saturation at a given fraction corrected for end effects, −
S_{i} : Normalized phase saturation, −
S_{1} : Normalized water saturation at inlet (Y = 1), −
S_{eq} : Normalized water saturation at which capillary pressure is zero, −
S_{r} : Reference scaled saturation (obtained if no end effects present), −
: Normalized water saturation averaged over the core, −
u_{i} : Darcy phase velocity, m/s
v_{i} : Interstitial velocity, m/s
Y : Scaled distance from outlet, −
Greek
σ_{ow} : Interfacial tension, N/m
λ_{i} : Phase mobility, 1/(Pa s)
Δp_{i}: Phase pressure drop, Pa
Δp_{i,ref}: Pressure drop without end effects, Pa
Indices
eq: Zero capillary pressure condition
r : Reference (no end effects)
1 Introduction
Relative permeabilities are an essential input to simulation when modeling multiphase flow processes, including petroleum recovery, carbon storage, and liquid invasion in hydraulic fracturing (Jeong et al., 2021; Juanes et al., 2006). Relative permeabilities describe the reduction of mobility for a flowing phase in presence of other phases and are modeled as function of saturation, although the saturation path, wettability, hysteresis, stress, and temperature and viscous coupling can affect them (Andersen et al., 2020a; Anderson, 1987; Bourbiaux and Kalaydjian, 1990). Their measurement should hence be performed under conditions close to the expected process taking place in the reservoir (Sidiq et al., 2017).
Relative permeabilities are commonly measured experimentally using the steady state technique, unsteady state technique or centrifuge. The steady state technique will be the focus of this work: It consists of injecting two fluids in different flow fractions (which can include fractions corresponding to single phase injection) and measuring the average saturation and the pressure drop over the core at steady state, when production rates and pressure readings have stabilized. If sufficient time has passed to reach steady state and the saturations are uniform, Darcy’s law can be applied to calculate relative permeabilities directly (Richardson et al., 1952). The unsteady state method considers injection of one phase to displace another where use is made of the transient production and pressure data before steady state is reached to calculate relative permeability (Johnson et al., 1959). With centrifuge, rotation at high speed allows calculating relative permeability for the less mobile phase from transient production data (Hagoort, 1980). All these methods are complicated by the presence of capillary pressure. This will be discussed only for the steady state method. Relative permeability can also be predicted from theoretical bundle of tubes models or digital rock models describing flow in representative porous structures (Nguyen et al., 2006; Valavanides, 2018).
In an open space there is no confinement or curvature of the fluid–fluid interface, setting capillary pressure equal zero at the producing end face of a core or an open fracture in a reservoir (Leverett, 1941). Pressure continuity forces capillary pressure profiles to converge to zero at the core outlet. Since there is a unique relation between capillary pressure and saturation (when saturation is changed monotonously) the steady state saturation profile will converge to the saturation giving zero capillary pressure. This results in nonuniform saturation profiles with most deviation at the outlet (Richardson et al., 1952). The capillary end effects can impact the calculation of relative permeabilities (Osoba et al., 1951). End effects cause steady state average saturation and pressure drop across the core to differ from what their values would be without end effects. Increasing the ratio of advective to capillary forces suppresses the end effect and makes the saturations more uniform and in accordance with the injected flow fraction. End effects can result in and explain apparently ratedependent relative permeabilities and flooding behavior, as observed experimentally (Alizadeh et al., 2007; Henderson et al., 1998; Jeong et al., 2021; Odeh and Dotson, 1985; Osoba et al., 1951; Rapoport and Leas, 1953). Chen and Wood (2001) found rateinsensitive relative permeabilities when their insitu measurements indicated negligible fluid accumulation at the outlet. Zou et al. (2020) used insitu imaging of a core plug steady state saturation profile to calculate relative permeability, but required an independent measurement of the capillary pressure function. This was achieved both via experimental and pore scale imaging based techniques. Rapoport and Leas (1953) showed that by increasing the core length, injected fluid viscosity and injection rate the end effects were less significant. However, reaching sufficiently high rates that end effects are negligible is not always practical due to limitations on core integrity, flow rate capacity, pressure reading, flow regime, and unrepresentative mobilization of residual droplets.
By varying the flow rate at a given flow fraction it is possible to assess the importance of the end effects and correct the measurements to representative saturations and pressure drops. Gupta and Maloney (2016) found that plotting pressure drop against rate was linear and gave a constant pressure drop at zero rate associated with end effects. Removing this extra pressure drop at each rate gave corrected consistent relative permeabilities. Further, they linked the extent of the end effect profile (how much of the core it covered) to how large the added pressure drop was relative to the corrected pressure drop. Plotting average saturation against this ratio could be extrapolated linearly to when the end effect had zero extent. The saturation and relative permeabilities obtained from these linear extrapolations are corrected for end effects and based on where the lines intercept. They hence called this the intercept method. Their main assumption was that the end effect saturation profile was limited to within the core. No relation was made to physical system properties, flow conditions or saturation functions. Their assumptions were mainly justified by numerical and experimental examples, rather than theory. A review of the intercept method was given by Reed and Maas (2018).
Andersen et al. (2017a) derived explicit analytical solutions for water flooding with end effects assuming specific saturation function correlations. They derived the intercept method theoretically, expressed using that average saturation and pressure drop divided by rate behaved linearly with inverse rate. They found how the saturation and pressure line slopes were related to input parameters of saturation functions. A capillary number was found determining the linear behavior with a critical value of 1 determining when the linear behavior was valid. They also could predict behavior when the end effects were outside the core. The model was validated experimentally by Andersen et al. (2020b) to determine both water relative permeability and capillary pressure. Huang and Honarpour (1998) also considered the waterflooding case, but under different saturation function assumptions. Their solutions were implicit, but could be used to correct end point relative permeability. Andersen and Zhou (2020) considered coinjection tests under the constraint of linear saturation functions and derived a capillary number incorporating all system parameters including the saturation functions. Virnovsky et al. (1995) demonstrated how sensitivity in average saturation and phase pressure drops to injection rate were analytically related to saturation functions. Andersen et al. (2020b) and Santos et al. (2021) history matched relative permeability and capillary pressure from multirate tests. It has also been demonstrated that saturation functions can be determined from cocurrent spontaneous imbibition experiments by varying the viscosity ratio systematically (Andersen, 2021; Andersen et al., 2019).
In this work we present mathematical relations describing end effects during steady state tests from fundamental assumptions and derive general conclusions regarding saturation profiles, parameters that affect saturation profile shape or just compress the profile, the intercept method, scaling, and dimensionless numbers.
Calculating effective relative permeabilities with Darcy’s law from data with end effects will be inaccurate, but combining such data from different total rates can allow accurate prediction. We derive a more intuitive and more general intercept method compared to the work of Gupta and Maloney (2016) based on established assumptions in core scale multiphase flow simulation. We show, for given fluids and injected flow fraction but varied total rate, that plotting steady state average saturation and inverse effective relative permeability against inverse total rate yields straight lines (at high total rates). The lines result in correct relative permeability and saturation points for that flow fraction when extrapolated to zero inverse total rate (infinite total rate). It is shown how the slopes of these lines are related to the saturation functions, allowing us to scale all data universally to the same straight lines. Two dimensionless numbers are obtained that directly express when end effects become important and how much in terms of how much (a) the average saturation is changed and (b) the effective relative permeabilities are reduced compared to values unaffected by end effects. A third dimensionless number is derived with a critical value of ~1 as a universal criterion for when the linear trends with inverse total rate (the intercept method) begin. The critical value is directly reflecting that the end effect profile has reached the inlet end of the core. When the profile extends beyond the inlet, the trends become nonlinear. All the dimensionless numbers are divided into a part depending on saturation functions, injected flow fraction and viscosity ratio and a second part containing constant known fluid and rock parameters such as core length, porosity, interfacial tension, etc. Methodologies for deriving relative permeability and capillary pressure consistently, even based on combining data from tests with different fluid and core properties are presented and demonstrated on two datasets from the literature.
2. Theory
2.1 Transport equations
The mathematical description of 1D incompressible and immiscible flow of oil (o) and water (w) in a porous homogeneous medium under negligible influence of gravity is given by mass balance and Darcy’s law, respectively:(1)
ϕ is porosity, s_{i} saturation of phase i = o, w, u_{i} Darcy velocity, K absolute permeability, λ_{i} mobility, k_{ri} relative permeability, μ_{i} viscosity, and p_{i} pressure. The saturations are dependent due to volume conservation, and the pressures are related by the capillary pressure function:(3)
The total Darcy velocity u_{T} is defined as:(4)
Also, λ_{T} is the total mobility, given by:(5)
It follows from adding the transport equations in (1) that:(6)
The water phase equation can then be expressed with variables u_{T}, s_{w} as:(7)
where f_{w} is the fractional flow function defined by:(9)
2.2 Boundary and initial conditions
Water and oil are injected simultaneously at the inlet x = 0 with a water flow fraction F (the water fraction of the total injected flux) and a total Darcy flux u_{T} (Fig. 1):(10)
Fig. 1 Illustration of the system including flow and boundary conditions, a typical end effect region and the relevant saturation interval at steady state. denotes the saturation where capillary pressure is zero, while denotes the saturation where the flow function f_{w} equals the injection flow fraction F. 
The injected water flux given u_{T} and F is then:(11)
The water flux (and that of oil) is composed of both an advective and capillary component, see (8). Hence, F does not correspond to f_{w} unless the capillary pressure gradient can be ignored. From (8) we write this boundary condition as:(12)
The outlet boundary condition is described by a zero capillary pressure (Leverett, 1941), which corresponds to a fixed outlet water saturation:(13)
where by definition P_{c} = 0.
2.3 Steady state
At steady state we have no changes with time in the system, i.e.:(14)
The phases are nonuniformly distributed due to the balance between advective and capillary forces. Given that time is not influential at steady state; in the following, water saturation and water pressure will be taken as functions of spatial coordinate alone: s_{w} = s_{w}(x) and p_{w} = p_{w}(x). Equation (7) can then be written as:(15)
At steady state the fluxes are uniform, i.e. the same amount of water and oil passes through every cross section, however the saturations and velocities can differ. Setting the water flux uniformly equal to that at the inlet, see (12), gives:(16)
Using that , we can solve (16) with respect to the saturation gradient:(17)
The water saturation gradient is thus dependent on the two phase mobilities, the capillary pressure curve, the injected water flow fraction F and the injection flux u_{T}. We can further introduce the interstitial total velocity, v_{T}, and dimensionless Leverett Jfunction (Dullien, 2012):(18)
The above equation can be integrated to find the saturation distribution starting from s_{w}(x = L) = . The saturation gradient will be nonzero until a saturation is reached such that:(20)
After which the saturation remains stable at . This is the state corresponding to negligible end effects. is found by solving (20). The pressure gradients of oil and water at steady state follow from (2) combined with (16):(21)
The above corresponds to Darcy’s law, where the water and oil fluxes are constant equal to u_{T}F and u_{T}(1 − F) and the mobilities vary along the core according to the steady state saturation distribution found from (19).
2.4 Scaled saturation profile
2.4.1 Derivation
Assume a domain where every saturation has a unique position. Equation (19) can then be solved by separation into a space coordinate integral and a saturation integral:(22)
Although the former integral is trivial, the latter in most cases requires numerical methods. Note that the latter saturation integral above has been expressed using normalized saturation S which provides the following relations:(23)
The parameters s_{wr}, s_{or} denote the critical saturations of water and oil, respectively, where their respective relative permeability is zero. Δs_{w} denotes the magnitude of the mobile saturation interval. S_{eq} is the normalized saturation where capillary pressure is zero.
Note that the terms f_{w} and f_{w}λ_{o} can both be written as functions of viscosity ratio, and that the latter term is inversely proportional to the geometric mean of viscosities:(26)
From this it is convenient to introduce the notations:(27)
N_{0} is a dimensionless capillary number (ratio of viscous to capillary forces) containing static or single phase flow parameters. This number is fixed if the same fluids, core and total rate are considered. N_{0} > 0 since v_{T} > 0 (flow in xdirection towards outlet). This leads to the solution form of interest:(28)
which is valid for all S between S_{eq} and S_{r}. S_{r} is the normalized saturation obtained at the flow fraction F without end effects:(29)
2.4.2 Shape characteristics
Consider two saturations S_{1}, S_{2} with positions given by (28). Their relative position is given by:(30)
This indicates that only the saturation functions, viscosity ratio and injected fraction determine the shape. The parameters in N_{0} will affect how compressed the profile is, but not its shape.
2.4.3 Length and range of end effect region
Based on the previous notation (27) the saturation gradient can be written from (19) as:(31)
Assume first that two phases are injected simultaneously (0 < F < 1) which means S_{r} ∉ {0, 1}. Further, consider saturations far from the outlet Y ≫ 0. The coefficient will then obey 0 < α_{0} < α(S) < α_{1} for some finite limits α_{0}, α_{1}. As distance from outlet increases the saturation will approach S_{r} according to the following Taylor approximation:(32)
Y will increase when S → S_{r}. We have considered a saturation region sufficiently close to S_{r} that α and can be considered constant as if evaluated at S_{r}. Assume a saturation S^{*} in that region with a finite difference from S_{r}. We integrate (32) and find that the spatial distance between the two saturations is infinite:(33)
S_{r} is only reached at infinite distance, in other words, it takes an infinite distance for end effects to vanish. This was exemplified by Andersen and Zhou (2020) for the special case of linear saturation functions, but is now proved in general.
On the other hand, assume single phase injection as given by F ∈ {0, 1}, which also implies S_{r} ∈ {0, 1}. The term f_{w}λ_{o} is the saturation dependent part of what is referred to as the capillary diffusion coefficient. It has parabolic shape and equals zero (only) at the points S = {0, 1}. Assume scaled saturations S approaching S_{r} and evaluate the saturation gradient:(34)
As both nominator and denominator in the first limit approach zero, we have applied L’Hopital’s rule. There exist parameter choices where could be zero or nonzero at the end points; using Corey relative permeabilities with exponents greater than 1 or equal to 1, respectively, will give such behavior Andersen et al. (2017a, 2020b) showed that the use of Coreytype functions for both relative permeability and capillary pressure resulted in a finite length end effect. Huang and Honarpour (1998) on the other hand used CoreyBurdine equations and obtained an infinite length.
The above derivations and examples show that the end effect will always have infinite length when two fluids are injected simultaneously. During single phase injection, the end effect will have either infinite or finite length depending on the saturation function correlations. We will however show that regardless of injection conditions and saturation function correlations we can define a practical length which produces the same results as if there was a finite length, and thus the intercept method can be derived. The intercept method was originally derived based entirely on the assumption that the end effect did not exceed the length of the core (Gupta and Maloney, 2016). It was found from several numerical examples that although the position of the saturation S_{r} extends to infinity, key saturation integrals converged as the integral limit approached S_{r}.
We have demonstrated the validity of the continuous saturation profile (28) from S_{eq} to S_{r} where a unique relation exists between Y and S. However, if S_{r} is finally reached at a specific Y(S_{r}), as can happen for the single phase injection case, all Y > Y(S_{r}) will have a constant saturation S_{r} which is necessary to account for as it impacts calculations of average saturation and pressure drop.
2.5 Dimensionless numbers based on cumulative end effect profile
Consider for simplicity that the porous medium extends infinitely beyond the core length Y = 1 such that saturations between S_{eq} and S_{r} have specified positions Y(S). The area between the straight line S = S_{r} and the graph Y(S) represents the cumulative amount of phase trapped by end effects as measured in displaceable pore volumes, n_{dpv}:(35)
The absolute sign on dS in the outer integral is used to produce a positive value regardless of whether S_{eq} is larger or less than S_{r}. As an example, if S_{r} = 0, S_{eq} = 1 and Y(S) = 1 (all the saturations 0 < S < 1 are positioned at the inlet end of the core), the area is 1 meaning n_{dpv} = 1. Whether the end effects have impact on the system will depend on how far the saturation profile is deviated from S_{r} throughout the core. n_{dpv} unfortunately cannot distinguish whether the accumulated phase is located within the core or not, but is a strong indicator of impact if the end effects are known to be within the core. The average length Y_{av} of the end effect profile (in core lengths) is:(36)
If Y_{av} = 0 there are no end effects while if Y_{av} = 1 the end effects are sufficiently strong to have a profile with average distance equal the length of the core. The saturation profile is nonuniform. Especially, the profile is locked at Y = 0 at S_{eq} so the saturations closer to S_{r} have distances greater than Y_{av}. Assume the profile Y(S) is approximated by an n’th order polynomial of S, such that Y_{av} is preserved and (S_{eq}) = 0:(37)
For n = 1 the profile is linear while larger n give the profile more curvature and results in the end effect region to terminate at position Y_{av}(n + 1):(38)
The parameter n can be selected to fit typical profiles or given an assumed value, but will be assumed fixed.
n_{dpv} and Y_{cee} are dimensionless capillary numbers expressing in different forms the ratio of capillary to viscous forces. They are derived from physical considerations and uniquely combine any set of system parameters, including the saturation functions. Large values ≫1 indicate strong end effects, low values ≪1 indicate negligible end effects and values of magnitude ≈1 are expected to give a transition in behavior. Especially, when Y_{cee} ≈ 1 the end effect profile should exactly reach the inlet. The two numbers are related by:(39)
Whether the end effects are limited to within the core or not is essential to the derivation of the intercept method (Andersen and Zhou, 2020; Andersen et al., 2017a, 2020b; Gupta and Maloney, 2016), suggesting Y_{cee} = 1 to be a criterion for when the method is valid.
2.6 Average saturation
2.6.1 General definition
In the following it will be assumed that the saturation profile extends to infinity, i.e. either coinjection is considered or we have single phase injection with a proper combination of saturation functions. The following spatial integrals are independent of this assumption, but their conversion to saturation integrals are not.
The core average saturation follows from integrating the saturation along the core. The integral can be converted into a saturation integral evaluated from the scaled outlet saturation S^{eq} to the scaled saturation S_{1} at the inlet Y = 1:(40)
In the above we have applied:(41)
which follows from (31). The inlet saturation S_{1} is unknown and depends on the extent of the saturation profile. It is found by solving the equation Y(S_{1}) = 1, equivalently:(42)
As seen, S_{1} depends on N_{0}, S_{eq} and the saturation functions. Especially, the saturation at this or any other specific location will depend on both the general profile shape and how compressed it is (the magnitude of N_{0}).
2.6.2 Intercept method
Assume now the saturation profile is approximated by introducing a saturation S^{*} such that:(43)
In other words, at a given fraction (1 − ε) of the saturation interval between S_{eq} and S_{r}, at a saturation S^{*} sufficiently close to S_{r}: the (infinite) end effect saturation profile is approximated to a finite end effect saturation profile stopping at Y(S^{*}). At greater distances there are assumed no end effects, S = S_{r}. The closer S^{*} is to S_{r} the better is the approximation. Further, assume that the length of the end effect is within the core: Y(S^{*}) < 1. This distance will be called Y^{*}. There is then a corresponding distance 1 − Y^{*} at the inlet without end effects. The average saturation in the core under such circumstances can be described by:(46)
is the average saturation in the end effect region (0 < Y < Y^{*}) and the position Y^{*} is known. By considering the saturation profile between S_{eq} and S^{*} we find the average saturation in that interval :(48)
As seen the average saturation in the end effect region is only a function of the saturation functions, flow fraction and viscosity ratio and not the parameters in N_{0}. The average saturation in the core however varies linearly with Y^{*}, the fraction of the core covered by end effects. Gupta and Maloney (2016) assumed this and verified it by running numerical simulations, but did not make a formal proof.
From the definition of Y^{*} and in (47) and (48) we can write (46) as:(49)
showing that the average saturation is linear with the inverse capillary number where the slope is a saturation integral and the intercept is the saturation without end effects.
The combination of the saturation integrals in (47) and (48) to one makes the final integral in (49) less sensitive to the choice of S^{*} and we can in fact let S^{*} → S_{r}. We then get a more correct (and less subjective) slope since at higher capillary numbers more of the profile will be within the core and should be accounted for.
The average (absolute) saturation can then be expressed as linear with the inverse capillary number 1/N_{0} with a saturation term slope C_{s} and, by expanding the capillary number expression, proportional with inverse velocity:(50)
The intercept is the corrected saturation without end effects .
2.7 Pressure analyses
2.7.1 Gradients
The pressure gradients of oil and water can be expressed in terms of the scaled distance from the outlet:(52)
2.7.2 Pressure drop
We define the pressure drop of a phase as the pressure at the inlet minus that at the outlet. We obtain these parameters by integration of the pressure gradients, expressed either as integrals over the positions of the core or saturation integrals, where we make use of knowing the saturation profiles:(54)
Assume for comparison that there were no end effects such that the pressure gradients (52) and (53) are constant and evaluated at S_{r}:(56)
Integrating from Y = 0 to 1, the corresponding pressure drops are:(58)
From the definition of S_{r} we have that f_{w}(S_{r}) = F which can be expanded to:(60)
and equivalently expressed for the oil phase using 1 − f_{w}(S_{r}) = 1 − F, we have:(61)
The two equations are both related to the total mobility and can therefore be combined:(62)
We thus see that the phase dependent terms in (56)–(59) are equal which implies that the phases have identical pressure gradients and pressure drop over the core in absence of end effects:(63)
We then also obtain identical phase pressure gradients once considering positions adequately far from the end effect zone. The ratios of pressure drop with end effects to pressure drop without end effects are:(65)
It is seen that N_{0} controls the impact of end effects on pressure drop for each phase, with the exception that the geometric viscosity in N_{0} is replaced by the phase viscosity (since N_{0} ∝ μ_{m}).
Of most interest from the pressure measurements are the data corresponding to the relative permeabilities without end effects k_{ri}(S_{r}). Assume that an “effective relative permeability” is calculated based on the measured pressure drop and injection conditions by direct application of Darcy’s law. On the other hand, the “true” relative permeability without end effects k_{ri}(S_{r}) would be obtained if the pressure drop was Δp_{ref}. This is related as follows:(67)
The ratio of pressure drops is then directly related to the ratio of relative permeability estimates:(68)
This implies the effective relative permeability approaches the true relative permeability when Δp_{i} approaches Δp_{ref} at high rates (small end effects):(69)
If the pressure drop is higher with end effects this leads to underestimation of the relative permeability , and vice versa.
2.7.3 The challenge of missing phase pressure data and its solution
Under normal circumstances we do not have access to the pressure drop of both phases, but measure only one phase pressure drop across the core. Virnovsky et al. (1995) state that during imbibition injection the pressures at the outlet are continuous, however only the nonwetting phase is continuous with the surroundings at the inlet. The experimentally measured pressure drop is therefore that of the nonwetting phase. Special inlet designs allowed pressure continuity and measurement of both fluid pressures (Virnovsky et al., 1998), but they are not common. An important question is then whether the available pressure drop still gives meaningful data to calculate the other phase’s relative permeability. Assume therefore that for a given phase the effective relative permeability is calculated using the pressure drop of the other phase. Equation (67) is then modified to:(70)
Since Δp_{j} approaches Δp_{ref} when end effects are negligible the correct relative permeability is still obtained, although the effective estimates with end effects are different when using the pressure drop of the other phase.
2.7.4 Intercept method
Again, assume the saturation profile is divided into a region 0 < Y < Y^{*} with end effects and a region Y^{*} < Y < 1 without end effects where S = S_{r}, as described before, where S(Y^{*}) ≈ S_{r}. Particularly, the end effect region is within the core and there exists a region without end effects. By integrating the pressure gradients over each interval, the pressure drops of each phase are:(72)
Dividing the pressure drops Δp_{i,ref} with the reference pressure drops Δp_{i} (obtained with no end effects) and using the definition of Y^{*} from (47) we can express the result in the following form:(74)
Again we note that the updated combined saturation integrals are less sensitive to S^{*} and we can let S^{*} → S_{r}. For both phase pressure drops scaled by the reference pressure drop, the result is linear trends with the inverse capillary number corrected for the phase viscosity times a saturation integral slope termed C_{w} and C_{o} for the respective phases. Expanding the capillary number shows that the relation is linear with inverse total velocity:(76)
The saturation term slopes are defined as:(78)
Equivalent to (68) we can use (76) and (77) to express the relations between effective and corrected relative permeabilities:(80)
which shows that if we plot the inverse of effective relative permeability against inverse capillary number or inverse velocity, we get a straight line trend with intercept at the ‘true’ relative permeability of each phase.
When calculating effective relative permeabilities for a given phase based on its own pressure drop data, we have shown that the intercept method follows where the inverse effective relative permeabilities plotted against inverse rate converge linearly to the correct inverse relative permeability without end effects. Assume now therefore that for a given phase the effective relative permeability is calculated using the pressure drop of the opposite phase. As shown in (71) this is directly related to that phase’s normalized pressure drop and combined with (76) and (77) we obtain:(82)
We have thus demonstrated that the intercept method still holds: the inverse relative permeability is linear with the inverse capillary number and converges to the correct inverse relative permeability. However, the slope of the linear data will correspond to the phase with continuous pressure.
2.8 Single phase injection
An important special case is single phase injection (fractions F equal 0 or 1) which is used to determine critical saturations and relative permeability end points. As before, average saturation is determined by (50). The relative permeability of the phase not flowing cannot be determined directly, but will by definition be zero at steady state for the corrected saturation . We now measure the pressure of the injected phase and its relative permeability is estimated using (80) if water is injected and (81) if oil is injected.
3 Results and discussion
Average (normalized) saturations can be calculated for any flow conditions based on (40) where the inlet saturation S_{1} is determined from (42). Absolute average saturation then follows using (23). Pressure drops are calculated from (54), (55) and reference pressure drop from (64).
For mixed or unknown wettability, at a given flow fraction, the lesswetting phase can be identified as the phase whose saturation increases upon increased injection total rate. The measured pressure drop at that fraction is then associated with the less wetting phase. From our model we therefore choose to only report the nonwetting phase pressure drop.
Effective relative permeabilities are calculated based on the nonwetting phase pressure drop. These data are valid under any flow conditions and can be compared directly with steady state solutions from any commercial simulator. In our case the nonwetting phase is oil and we use (67) for effective oil relative permeability and (70) for effective water relative permeability . An important exception is water flooding where the capillary pressure becomes zero or negative along the core and water will have higher pressure drop than the oil.
In addition, the analytical solutions corresponding to the intercept method are presented. They are valid only at high rates and will expectedly differ from the general solution at low rates. Selected effective relative permeability points are also calculated based on Sendra v2018.2.5 for validation.
3.1 Input data
Reference case input parameters are displayed in Table 1. As indicated, the oilwater flow parameters were based on measurements from Kleppe and Morse (1974) from a Berea sandstone core. Correlations with sufficient parameters to match the data were selected, such as the Andersen et al. (2017b) Jcorrelation for scaled capillary pressure and an extended Corey correlation (Brooks and Corey, 1964) for relative permeability with Corey exponents linearly dependent on saturation:(84)
Reference input parameters in the simulations. The parameters are based on Kleppe and Morse (1974, January). To scale the capillary pressure, interfacial tension was assumed to be 21 mN/m. The core length was assumed to be 10 cm.
All the functions were expressed using the normalized saturation S from (23). The Jfunction derivative, which is central in the solution is easily evaluated explicitly as:(87)
The relative permeabilities and scaled and unscaled capillary pressure can be seen in Figure 2 comparing correlations with the experimental data. Fractional flow functions with three choices of oil viscosity are shown in Figure 3.
Fig. 2 Input saturation functions based on Kleppe and Morse (1974, January): oil and water relative permeability (left) and scaled and absolute capillary pressure (right). The experimental points are shown together with the correlations based on (84)–(86). 
Fig. 3 Fractional flow functions f_{w} for different oil viscosities (the base value of 2.3 cP and increased values 5 and 25 times higher) against water saturation in lin–lin plot (left) and semilog plot (right). 
3.2 Validation
The analytical solution was validated by comparison of results with full transient numerical simulations from the commercial software Sendra v2018.2.5. This program solves the flow equations fully implicit. The core was discretized using 500 grid cells to capture steep saturation gradients. At each injection condition the program was run until a visibly flat production profile was observed, indicating steady state.
Using the base case parameters, two flow fractions are considered: F = 0.50 and F = 0.05 and injection rate is varied for both cases using 0.1, 1, 10 and 100 pore volumes per day (PV/d). The resulting steady state saturation profiles are shown in Figure 4 comparing the analytical and numerical solutions. For comparison, the saturation integrals in the derived methods were calculated using 10 000 saturations in the interval between S_{r} and S_{eq}, sufficient for the integrals to converge.
Fig. 4 Validation of the analytical model by comparing results from the analytical solution and a commercial software (Sendra) at high (a) and low (b) injected fraction (F = 0.5 and 0.05, respectively) at different total injection rates, as indicated in the legend, measured in pore volumes per day. 
For all rates and flow fractions the saturation distribution solutions from the analytical solution (28) and the commercial software (Sendra) overlap, verifying our procedure. All the saturation distributions approach the saturation at Y = 0. For a given flow fraction the saturation distributions approach the saturation such that f_{w} = F. Considering Figure 3 we see that when F = 0.5 and F = 0.05, then and ≈ 0.475, respectively, which is consistent with the profiles in Figure 4. For high rates the distributions stabilize at within the core a certain distance from the outlet. For low rates (such as 0.1 PV/d in these examples), however, the capillary forces dominate and shift the saturations towards and no saturations in the core are close to . Direct application of Darcy’s law to calculate relative permeability at steady state assumes which we see is only reasonable at very high rates.
3.3 Profile shapes
The analytical solution predicts that variation of N_{0} only compresses or expands the saturation profile as long as the saturation functions, viscosity ratio and flow fraction are not changed, see (44). This is demonstrated by varying the core length L, mean viscosity μ_{m}, total velocity v_{T} and permeability K (which all are directly involved in N_{0}, see (27)) around the base case with F = 0.05 and rate 1 PV/d to give same change in N_{0} (an increase of N_{0} by a factor 5). The corresponding saturation profiles are shown in Figure 5a together with the reference case. As expected, all the parameter variations increasing N_{0} by the factor 5 result in the same scaled profile which is only a compression of the base profile. For example, the saturation s_{w} = 0.5 has position Y ≈ 0.9 in the base case and = 0.18 in the new cases.
Fig. 5 Saturation profiles Y(s_{w}) calculated using (24) which consists of a dimensionless number N_{0} and a saturation integral. Variation of parameters in N_{0} (a) that increase N_{0} by the same factor from the base case, but do not change the saturation integral, result in an identical new case corresponding to a compression of the base case profile. Changing viscosity (ratio) or flow fraction (b) affect the saturation integral and can change the saturation range and profile shape in general. Compressing these profiles would not cause them to overlap. 
In comparison, in Figure 5b we see that variations in individual viscosities (thus changing the viscosity ratio) and flow fraction affect both magnitude and shape of the profiles as they do not overlap upon scaling. In particular, different values of the end point are obtained in each case.
3.4 Determining the nparameter for approximate profiles
The n parameter appears in the scaling number Y_{cee}. To find a good choice for n we want to see that the actual numerical profiles are well approximated by the polynomial function. This was considered by varying the injected fraction F (with values 0.01, 0.05 and 0.5) and the oil viscosity μ_{o} (with values 1, 5 and 25 times the reference value of 2.3 cP) to get distinct profile shapes. The normalized saturation profiles S(Y) (blue curves) were plotted together with the corresponding approximated profiles (orange curves) in Figure 6 for n = 5 which was found to be a good choice. The profiles are plotted over the range Y = 0 to 10Y_{av}. The positions Y_{av} and Y_{cee} = (n + 1)Y_{av} are marked with vertical dashed lines (left and right, respectively). Also a dashed/dotted line indicating the normalized saturation S_{r} (obtained at infinite distance) is shown. The main observations are:
The major part of the saturation profile is located behind Y_{cee} including 90–95% of the saturations closest to S_{eq}.
Saturations located at Y > Y_{cee} are very close to S_{r}.
We can state end effects to be severe when Y_{cee} > 1 since then all saturations deviate from S_{r} more than 5% the magnitude of the end effect saturation interval S_{r} − S_{eq}.
Fig. 6 Normalized saturation profiles from the analytical solution (in blue) for different injected flow fractions and oil viscosities and their corresponding approximated saturation profiles assuming polynomial shape for n = 5. The profiles are plotted against Y over the range 0 to 10Y_{av} where Y_{av} and Y_{cee} = (n + 1)Y_{av} are marked with dashed vertical lines left and right, respectively. 
In terms of the impact the fraction and viscosities have, we see that when the fraction F is lowered, S_{r} decreases farther from S_{eq}. This expands the saturation range and the extent of end effects. We also see that increasing the oil viscosity alters the mobility ratio to reduce S_{r} and expands the saturation range. Although this mechanism increases the end effect, the increased viscous forces work to compress the profile. At low F the increase in oil viscosity has more impact on compressing the end effect profile.
3.5 Intercept method
3.5.1 Demonstration of analytical solution/intercept method
Consider the reference case for the four combinations of flow fractions F = 0.05 and 0.01 and oil viscosities 2.3 cP and 5 × 2.3 cP. Average saturation and inverse effective oil and water relative permeabilities, and , are plotted against inverse velocity (Fig. 7). Flow rates are used in the range 5 PV/d to 1000 PV/d. The following is seen:
For all cases, once the rates are sufficiently high, the results overlap completely with the linear analytical solution. As , the solutions reach the correct saturation and relative permeability points without end effects, as marked on the yaxis with large symbols.
At low rates the two solutions do not overlap. The analytical solution is only valid when the end effect profile has not reached the outlet significantly which is less likely at low rates.
When plotting the relative permeability data on two scales such that two points overlap, all the points overlap, i.e. the profile shapes are identical.
Since the trends are linear, and the linear trends appear valid over a wide range of rates, only a few rates are necessary to determine the corrected values.
The above statements were validated for average saturation and relative permeability by simulations with the software Sendra for the case F = 0.01, μ_{o} = 2.3 cP.
Fig. 7 Plots of average saturation vs. inverse total rate (left) and plots of inverse relative permeabilities of oil and water vs. inverse total rate (right) under different flow fractions and oil viscosities. The results are calculated both using the general solutions (50), (54) and (55) which are valid under all conditions, and using the analytical solution (AnSol) in (50) and (80)–(83). For the case F = 0.01, μ_{o} = 2.3 cP simulation results from the commercial (COM) software Sendra are added for comparison. 
3.5.2 Scaling of end effects with dimensionless numbers
We have seen that expressing average saturation and inverse relative permeabilities as function of inverse velocity results in linear relations. The velocity dependent terms directly express the impact of end effects on those properties and are reasonable to apply as dimensionless numbers for scaling the behavior.
In Figure 8a we plot against the number based on different flow fractions (F = 0.9–0.005) and injection rates and in Figure 8c we do the same varying oil viscosities (from 5^{0} to 5^{4} times the reference value) and injection rates. For a given flow fraction, C_{s} is constant while N_{0} varies with rate. Oil viscosity changes both terms.
Fig. 8 Comparison of the general solutions (full or dash/dotted lines) to the analytical solution (straight dashed line) in terms of in (a) and (c) and in (b) and (d) for changes in flow fraction F in (a) and (b) and changes in oil viscosity in (c) and (d). All solutions fall on the same straight line at low dimensionless number (high rates) which scale the data. 
At high rates is identical to , which follows directly from the analytical solution, see (50), and all cases fall on this straight line. The dimensionless number thus scales saturation end effects and accounts for all relevant input parameters. For low rates the cases deviate from the straight line with different trends. The number directly expresses the change in average saturation due to end effects. As end effects also can reduce the average saturation, the straight line can continue to negative values. Values of can be considered ‘large’ if it changes the average saturation by 0.01 or more (for comparison the mobile saturation range in our data is 0.31) and negligible for .
We also plot against in Figure 8b for different flow fractions and Figure 8d for different oil viscosities. From (80) to (83) we see that the two expressions are identical and form a straight line from zero at high rates (low ). At low rates the data deviate from the straight line.
On the straight line, the term expresses how much the effective relative permeability is reduced compared to the corrected relative permeability as given by . If , they are identical, while if nonzero they differ. A value of greater than 1 can be considered large since then . Since the measured pressure drop is in the less wetting phase which has the highest pressure drop, will be positive and the effective relative permeabilities will be lower than the corrected relative permeabilities, i.e. . This is seen in the graph as .
If conditions were such that water was nonwetting at a given fraction (in a mixedwet or oilwet system) then we should modify the plot to show against .
An interesting observation is that when the curves deviate from the straight line, they all fall below it, indicating that the impact is less than if the entire end effect profile was within the core. That is reasonable since the profile that is not within the core does not contribute to change the observations, and since at lower rates more of the profile is outside the core, the added impact on the saturation and relative permeability becomes less.
Note also that the dimensionless number that scales saturation, , is different from the one scaling relative permeability, or . As mentioned, the dimensionless number which is used for scaling relative permeability is also determined by which phase is more wetting at the given fraction (i.e. whether is above or below which further depends on wettability, flow fraction and viscosity ratio).
To see more systematically when the solutions deviate from the analytical solution (the linear model) we plot and for the aforementioned cases in Figure 9, using the same dimensionless numbers and , respectively. When the analytical solution is valid, these terms are 0 and 1, respectively. This is the case at low values of the respective dimensionless numbers. The onset of deviation in saturation becomes visible for 0.005 < < 0.02 for the cases with variation in F and 0.008 < < 0.03 for the cases with variation in oil viscosity. The onset of deviation in relative permeability trends from 1 is visible for 0.8 < < 5 in the cases with variation in F and 0.8 < < 3 for the cases with variation in oil viscosity. The onset thus occurs at values of the dimensionless numbers varying by an order of magnitude. Further, both dimensionless numbers appear directly linearly in their separate relations (for saturation or relative permeability). However, their magnitudes at which we observe deviation from the linear relations are very different, spanning three orders of magnitude. Hence, although the dimensionless numbers are excellent at scaling the end effects, they do not accurately determine when the linear analytical solutions cease to be valid.
Fig. 9 The difference between general solutions and the analytical solution in terms of plotted against C_{s}/N_{0} in (a) and (c) and plotted against in (b) and (d) for changes in flow fraction F in (a) and (b) and changes in oil viscosity in (c) and (d). The general solution begins to deviate from the analytical solution at a range of dimensionless numbers. 
However, we found that Y_{cee} is a good measure of when significant changes in saturation have occurred along the entire core and previously found an appropriate value of n = 5 to apply in this number from matching saturation profiles, see Figure 6. We now present the above results of and plotted against Y_{cee} in Figure 10. For all cases and both the saturation and relative permeability data we see that the onset of deviation starts at Y_{cee} ≈ 1 with very little variation (a factor 1.2, i.e. 0.08 orders of magnitude). As we have reduced the span from 1 order of magnitude for individual terms and from 3 orders of magnitude considering both terms and to 0.08 orders we can claim that Y_{cee} < 1 is a reliable criterion for when the intercept method is valid.
Fig. 10 The difference between general solutions and the analytical solution in terms of in (a) and (c) and in (b) and (d) for changes in flow fraction F in (a) and (b) and changes in oil viscosity in (c) and (d). All results are plotted against Y_{cee} and demonstrate that the analytical solution is valid for Y_{cee} < 1. 
3.6 Appearance of ratedependent relative permeabilities
End effects can affect both the steady state average saturation and effective relative permeability points, which deviate from the true relative permeability functions. Effectively this can result in an appearance of ratedependent relative permeability functions meaning that if they are measured at a fixed rate, the rate will affect which curves are obtained and not offer information about their correctness.
Based on the reference case relative permeability functions and capillary pressure as input, effective relative permeabilities were calculated and plotted against average saturation. For a given injection total rate, 40 evenly spaced fractions were used between F = 0.005 and 0.995, each producing a relative permeability point. This was performed for nine injection total rates ranging from 0.1 to 1000 PV/d varied by the same factor 10^{0.5}. Since not only rate, but the dimensionless number N_{0} determines the impact of end effects, the effective relative permeability curves were labeled by the corresponding value of N_{0} and are shown in Figure 11. High rate corresponds to high N_{0} and opposite. The input relative permeabilities, which are the true curves of the system, are included for comparison with the effective curves.
Fig. 11 Effective relative permeabilities and saturations calculated for fractions F between 0.005 and 0.995 for injection total rates varying from 0.1 to 1000 PV/d. The calculated number N_{0} is shown, where high N_{0} corresponds to high injection rate. Each curve is based on a fixed N_{0}, and each saturation on a given curve is obtained from one value of flow fraction. The reference (“true”) input relative permeabilities and capillary pressure curves were used and were only approximated by the effective curves at sufficiently high N_{0}. End effects are important at low N_{0} and shift the saturations towards = 0.61 and reduce the effective relative permeabilities. 
Notably, at low N_{0} (low total rates) the average saturations are shifted towards = 0.61, leading to narrower saturation ranges. Each calculated relative permeability point is reduced by the high pressure drop in the nonwetting phase. At high N_{0} (high total rates) the relative permeability curves converge and become independent of injection total rate (and N_{0}). They then also agree with the input curves (“true rel perms”) used to generate the data.
We do not argue that relative permeabilities physically cannot exhibit ratedependence. For example, micromodel observations indicate the possible existence of a flow regime dependent of the conventional capillary number (velocity times viscosity divided by interfacial tension) (Valavanides, 2018).
3.7 Calculation of relative permeability and capillary pressure from experiments
We consider steady state experiment data where two phases are injected at different flow fractions F^{j} and at each fraction different total velocities are applied. Potentially such tests can be applied using different cores and even fluids, under the assumption that the saturation functions remain the same. For example, different core tests can be implemented at constant total velocity using the same schedule of injected fractions (Henderson et al., 1998; Virnovsky et al., 1998). The cores used can have different length L, porosity ϕ, permeability K. Use of different fluids (or just measuring a different value for the same fluids in a different test) can change viscosities μ_{i} and interfacial tension σ_{ow}. Average saturation and pressure drop Δp^{j,k} are measured at each fraction and velocity.
3.7.1 Simple procedure: Calculation of only relative permeability with constant fluid and rock properties
Convert the pressure drop data to effective relative permeability measurements .
For each fraction F^{j} plot average saturation against inverse total rate and inverse effective relative permeabilities against inverse total rate . Determine the lines through the data:
Report the intercepts . The slopes are not used.
The pairs for i = o, w are accurate relative permeability points. Match suitable relative permeability correlations (and critical saturations) to fit them.
3.7.2 Advanced procedure: Calculation of both relative permeability and capillary pressure from tests with distinct properties
For each fraction F^{j} and total velocity , calculate the number and the viscosity ratios . The involved parameters are standard.
For each fraction F^{j} plot average saturation points against . Determine the line through the data (all or those at lowest ) by regression and report the slope and intercept as:
For each fraction F^{j} determine which phase is more wetting. If average water saturation decreases at lower , water is wetting at that fraction. If opposite, oil is wetting. The measured pressure drop corresponds to the less wetting phase.
Convert the pressure drop data to effective relative permeability measurements :
Assume first the fluids have constant viscosities. For each fraction F^{j} and rate plot inverse effective relative permeabilities against . Determine the line through the points (all or those at lowest ) and report the slope and intercept as:
Both relative permeability lines will have the same slope . Depending on which phase is wetting at the given fraction we determine either or :(94)
If the viscosities vary between tests, the data must be plotted against for the fractions where water is less wetting and against where oil is less wetting.
The pairs for i = o, w are accurately determined relative permeability points.
Match suitable relative permeability correlations (and critical saturations) to the relative permeability points.
The slopes indicate how strongly capillary forces affect saturation and pressure drop at the relevant saturation ranges to and constrain the choice of capillary pressure functions. Their values only depend on the saturation functions and are independent of fluid and core properties.
Assume a capillary pressure correlation with tuning parameters. Based on the calculated relative permeability correlations, select tuning parameters and calculate slope constants . The tuning parameters giving least relative error between predicted and measured values result in the capillary pressure Jfunction best explaining the data.
3.7.3 Interpretation of experimental data I
Virnovsky et al. (1998) conducted three oilwater steady state tests on a Berea sandstone core. The core was initially fully water saturated and oil fractions were increased. A fixed total rate was used in each test of 0.2, 0.5 or 5 cc/min. The absolute permeability changed between the three tests, hence we used N_{0} to account for the differences. From the three velocities at each fraction we follow the ‘advanced’ procedure. Average saturations and inverse effective relative permeabilities were plotted against in Figures 12a–12d, respectively to find the intercepts and the slope constants C_{s}, C_{o}, C_{w} which are listed in Table 2. Linear trends in the points were clearly visible for most cases.
Fig. 12 Interpretation of experimental data (crosses) from Virnovsky et al. (1998) where average saturation (a) and inverse effective relative permeabilities in (b) and (c) and in (d) and (e) are plotted against . Straight lines are drawn through the data for each flow fraction F to find the intercept (marked with circle points) corresponding to saturations and relative permeability points corrected for end effects. The slopes of the lines provide values of C_{s}, C_{o}, C_{w} at each fraction. 
Line analysis of data from Virnovsky et al. (1998) giving corrected saturations and relative permeability points from the intercepts, and slope values C_{s}, C_{o}, C_{w} to be used for derivation of capillary pressure for each fraction F.
For the 5 highest fractions F water saturation increased with rate, indicating capillary forces trapped oil then, while for the lowest 3 fractions water saturation decreased with increased rate, indicated capillary forces then trapped water. Equivalently C_{s}, the slope of average saturation vs. , was negative for the highest fractions and positive for the lowest fractions. Although there was some uncertainty in the sign for some central fractions it seems clear that was located in the intermediate saturation range.
As described in the theory, we can scale end effects by plotting against and against (when oil is less wetting) and against (when water is less wetting). This was done for all the data points and is shown in Figure 13. Average saturation, water relative permeability and oil relative permeability fall on a straight line in their respective plot.
Fig. 13 Scaling of the data from Virnovsky et al. (1998) by plotting against , (b) and against (when oil is less wetting) and against (when water is less wetting). 
The relative permeability points were fitted with the correlation (85) and (86), see Figure 14a. Next, the Jfunction was determined, see Figure 14b, using the correlation (84) by optimizing the match of slope parameters: The C_{s} values were well matched, see Figure 15a, while it was more difficult to match the C_{o}, C_{w} parameters which at optimum were roughly 2–10 times lower than those derived from the measurements, see Figure 15b. This means that the model with the tuned curves would capture well what the average saturation and effective relative permeability values would approach at each fraction, and that the slope of average saturation with inverse total rate would vary as observed experimentally. However, the effective relative permeabilities would change less in the model than observed experimentally.
Fig. 14 Relative permeabilities in (a) with points based on intercepts from the steady state measurements from Virnovsky et al. (1998) and the correlations (full lines) using (85) and (86) that best fit the points. The best fitting capillary pressure correlation is shown in (b) as scaled Jfunction and in mbar based on matching the slope parameters from the steady state measurements. The curve is compared to a primary drainage curve measured with centrifuge on the core used in Virnovsky et al. (1998). Our derived Jfunction is more consistent with the experimental data as reflected in a higher residual water saturation and the presence of negative capillary pressures. 
Fig. 15 Optimal match of the slope coefficients C_{s}, C_{w}, C_{o} after tuning the Jfunction. 
Near identical Jfunctions resulted from different optimizations, indicating that the data and approach give consistent results. The optimal curve parameters and the system parameters are found in Table 3. Comparison with a measured drainage capillary pressure curve indicated similar magnitude and shape, however negative capillary pressures at high saturations resulted from the model, required to explain the trends in saturation with rate. Also, a higher residual water saturation was found which seems more consistent with the flooding data.
System parameters from Virnovsky et al. (1998); Relative permeability correlation parameters to fit the corrected relative permeability data; Jfunction correlation parameters to fit the slope parameters from the same dataset. Experiments marked A, B, C with different velocities reported different permeability.
The model can under any conditions only be sensitive to saturations on the interval between and . As was located centrally, the extreme values corresponding to the lowest and highest fractions gave the total saturation range where we have accurate information. As seen in Table 2 and Figure 14, this is between saturations 0.326 and 0.967 which is practically the entire mobile saturation range.
3.7.4 Interpretation of experimental data II
Henderson et al. (1998) measured relative permeability to gas and condensate at 7 flow fractions between 0.048 and 0.29 (corresponding to condensategasratio 0.05–0.4) and 4 total velocities on Berea sandstone with presence of connate water saturation. Condensate c and gas g were treated as ‘oil’ and ‘water’ in the above theory. Plotting average condensate saturations and inverse effective relative permeabilities and against in Figure 16 indicates that the points at the three highest rates (lowest ) aligned well on straight lines. This was used to obtain the intercepts and resulting points for i = c, g, as well as the slope parameters C_{s}, C_{c}, C_{g}. All these values are listed in Table 4 for their respective flow fractions.
Fig. 16 Interpretation of experimental data (crosses) from Henderson et al. (1998) where average condensate saturation (a) and inverse effective relative permeabilities in (b) and in (c) are plotted against . Straight lines are drawn through the high rate data for each flow fraction F to find the intercept (marked with circle points) corresponding to saturations and relative permeability points corrected for end effects. The slopes of the lines provide values of C_{s}·C_{o}, C_{w} at each fraction. 
Line analysis of data from Henderson et al. (1998) giving corrected saturations and relative permeability points from the intercepts, and slope values C_{s}, C_{c}, C_{g} to be used for derivation of capillary pressure for each fraction F. Condensate was less wetting than gas for all fractions, hence C_{g} was not reported for any fractions.
For all the data, the points at the lowest rate (highest ) that did not fall on the line, had less impact on the saturation or relative permeability than if they had been on the line (the trends flattened). This is consistent with what was shown in Figure 8 where the end effects’ added impact lessens as more of the profile is outside the core.
Also, at each fraction, increased velocity at a fixed fraction increased the condensate saturation, or equivalently condensate saturation decreased with higher , quantified by negative C_{s}. This indicated that capillary forces trapped gas under all conditions and gas wetted the rock more than condensate for the considered fractions.
Scaling the saturation and relative permeability points is shown in Figure 17 where we again see that all the high rate points fall well on the same straight line when plotted against the relevant dimensionless numbers. The low rate results fall between the line and the horizontal axis indicating that the impact is less than if they had continued on the line.
Fig. 17 Scaling of the data from Henderson et al. (1998) by plotting against , (b) and against (when gas is less wetting) and against (when condensate is less wetting). 
The intercept saturation and relative permeability points were collected and matched with correlations, as shown in Figure 18a. The slope parameters were then matched by tuning the Jfunction, seen in Figure 18b. Only the saturation range where the curves are reliable is shown. That includes the range between the lowest and highest (0.397–0.507), but also down to = 0.344 since all saturation profiles end at that saturation. The optimal match between the slope parameters from the tuned model and the values from the observations is shown in Figure 19. Both saturation and relative permeability slope parameters are well matched, with almost all values from the tuned model being less than 30% different from the values from the measurements. System parameters and correlation parameters from matching the data are listed in Table 5.
Fig. 18 Relative permeabilities in (a) with points based on intercepts from the steady state measurements from Henderson et al. (1998) and the correlations (full lines) using (85) and (86) that best fit the points. The best fitting capillary pressure correlation is shown in (b) as scaled Jfunction and in mbar based on matching the slope parameters from the steady state measurements. 
Fig. 19 Optimal match of the slope coefficients C_{s} (a) and C_{c} (b) after tuning the Jfunction. C_{g} was not relevant. 
System parameters from Henderson et al. (1998); Relative permeability correlation parameters to fit the corrected relative permeability data; Jfunction correlation parameters to fit the slope parameters from the same dataset.
Henderson et al. (1998) interpreted their data as capillary number dependent relative permeabilities rather than being affected by capillary end effects. They suggested that the pressure drop across the core was orders of magnitude higher than the capillary pressure that could cause end effects. Although we cannot conclude that end effects were present or not, they explain the observations very consistently in that the data fall on the predicted straight lines at high rates, fall on the expected side of the line at low rates and we can obtain saturation functions that consistently explain all the data. Gupta and Maloney (2016) interpreted the same dataset to find corrected relative permeabilities.
4 Conclusion
A theory is derived from fundamental assumptions describing how capillary end effects affect average saturation and relative permeability calculation in steady state experiments. The method is valid for all saturations functions and wetting states.
We derive an intuitive and general ‘intercept’ method for correcting steady state relative permeability measurements for capillary end effects: plotting average saturation and inverse relative permeability against inverse rate will give linear trends at high rates and result in corrected relative permeability points when extrapolated to zero inverse rate (infinite rate). This is a formal proof and generalization of the method proposed by Gupta and Maloney (2016).
It is shown how the slopes of the lines are related to the saturation functions allowing to scale all data for all conditions to the same straight lines. Two dimensionless numbers are obtained that directly express when end effects become important in terms of (a) how much the average saturation is changed: C_{s}/N_{0} and (b) how much effective relative permeabilities are reduced .
A third dimensionless number Y_{cee} (expressing the scaled front position of the end effect region) is derived directly stating that the end effect profile has reached the inlet of the core significantly when Y_{cee} = 1. This number acts as a universal criterion for when the linear (intercept method) behavior begins. The intercept method is valid when Y_{cee} < 1.
All the dimensionless numbers contain a part depending only on saturation functions, injected flow fraction and viscosity ratio and a second part N_{0} containing constant known fluid and rock parameters such as core length, porosity, interfacial tension, etc. The former parameters determine the saturation range and shape of the saturation profile, while the number N_{0} determines how much the profile is compressed towards the outlet.
End effects cause the average saturations to be shifted towards (the saturation at which capillary pressure is zero) and the calculated relative permeabilities from pressure drop and Darcy’s law to be reduced in magnitude compared to the true relative permeabilities. Since the shift of the effective curves depends on the viscous vs. capillary forces, using different rates will result in different effective curves when capillary end effects are significant. This gives a false impression of ratedependent relative permeabilities.
Methodologies for deriving relative permeability and capillary pressure systematically and consistently, even based on combining data from tests with different fluid and core properties, are presented and demonstrated on two datasets from the literature.
Even with access to only one phase pressure drop of two injected phases it is shown that the intercept method holds for both phase relative permeabilities. This is the standard regarding what pressure data is measured.
Acknowledgments
The author acknowledges the Research Council of Norway and the industry partners, ConocoPhillips Skandinavia AS, Aker BP ASA, Vår Energi AS, Equinor ASA, Neptune Energy Norge AS, Lundin Norway AS, Halliburton AS, Schlumberger Norge AS, and Wintershall DEA, of The National IOR Centre of Norway for support.
References
 Alizadeh A.H., Keshavarz A., Haghighi M. (2007) Flow rate effect on twophase relative permeability in Iranian carbonate rocks, in: SPE Middle East Oil and Gas Show and Conference, March 11–14, 2007, Manama, Bahrain. Society of Petroleum Engineers. [Google Scholar]
 Andersen P.Ø., Standnes D.C., Skjæveland S.M. (2017a) Waterflooding oilsaturated core samples – Analytical solutions for steadystate capillary end effects and correction of residual saturation, J. Pet. Sci. Eng. 157, 364–379. [Google Scholar]
 Andersen P.Ø., Skjæveland S.M., Standnes D.C. (2017b) A novel bounded capillary pressure correlation with application to both mixed and strongly wetted porous media, in: Abu Dhabi International Petroleum Exhibition & Conference, November 13–16, 2017, Abu Dhabi, UAE. Society of Petroleum Engineers. [Google Scholar]
 Andersen P.Ø., Brattekås B., Nødland O., Lohne A., Føyen T.L., Fernø M.A. (2019) Darcyscale simulation of boundarycondition effects during capillarydominated flow in highpermeability systems, SPE Reserv. Eval. Eng. 22, 02, 673–691. [Google Scholar]
 Andersen P.Ø., Nesvik E.K., Standnes D.C. (2020a) Analytical solutions for forced and spontaneous imbibition accounting for viscous coupling, J. Pet. Sci. Eng. 186, 106717. [Google Scholar]
 Andersen P.Ø., Walrond K., Nainggolan C.K., Pulido E.Y., Askarinezhad R. (2020b) Simulation interpretation of capillary pressure and relative permeability from laboratory waterflooding experiments in preferentially oilwet porous media, SPE Reserv. Eval. Eng. 23, 01, 230–246. [Google Scholar]
 Andersen P.Ø., Zhou Y. (2020) Steady state relative permeability experiments with capillary end effects: analytical solutions including derivation of the intercept method, J. Pet. Sci. Eng. 192, 1–9. [Google Scholar]
 Andersen P.Ø. (2021) Early and late time analytical solutions for cocurrent spontaneous imbibition and generalized scaling, SPE J. 26, 1, 220–240. [Google Scholar]
 Anderson W.G. (1987) Wettability literature survey part 5: The effects of wettability on relative permeability, J. Pet. Technol. 39, 11, 1453–1468. [Google Scholar]
 Bourbiaux B.J., Kalaydjian F.J. (1990) Experimental study of cocurrent and countercurrent flows in natural porous media, SPE Reserv. Eng. 5, 03, 361–368. [CrossRef] [Google Scholar]
 Brooks R.H., Corey A.T. (1964) Hydraulic properties of porous media. Hydrology Papers, No. 3, Colorado State U., Fort Collins, Colorado. [Google Scholar]
 Chen A.L., Wood A.C. (2001 September) Rate effects on wateroil relative permeability. In Proceedings of the International Symposium of the Society of Core Analysts, 17–19 September 2001, Edinburgh, Scotland, pp. 17–19. [Google Scholar]
 Dullien F.A. (2012) Porous media: Fluid transport and pore structure, Academic Press, San Diego, California. [Google Scholar]
 Gupta R., Maloney D.R. (2016) Intercept method – A novel technique to correct steadystate relative permeability data for capillary end effects, SPE Reserv. Eval. Eng. 19, 02, 316–330. [Google Scholar]
 Hagoort J. (1980) Oil recovery by gravity drainage, Soc. Pet. Eng. J. 20, 03, 139–150. [Google Scholar]
 Henderson G.D., Danesh A., Tehrani D.H., AlShaidi S., Peden J.M. (1998) Measurement and correlation of gas condensate relative permeability by the steadystate method, SPE Reserv. Eval. Eng. 1, 02, 134–140. [Google Scholar]
 Huang D.D., Honarpour M.M. (1998) Capillary end effects in coreflood calculations, J. Pet. Sci. Eng. 19, 1–2, 103–117. [Google Scholar]
 Jeong G.S., Ki S., Lee D.S., Jang I. (2021) Effect of the flow rate on the relative permeability curve in the CO_{2} and brine system for CO_{2} sequestration, Sustainability 13, 3, 1543. [Google Scholar]
 Johnson E.F., Bossler D.P., Bossler V.O. (1959) Calculation of relative permeability from displacement experiments, Trans. AIME 216, 01, 370–372. [Google Scholar]
 Juanes R., Spiteri E.J., Orr F.M. Jr, Blunt M.J. (2006) Impact of relative permeability hysteresis on geological CO_{2} storage, Water Resour. Res. 42, 12, 1–13. [Google Scholar]
 Kleppe J., Morse R.A. (1974, January) Oil production from fractured reservoirs by water displacement, in: Fall Meeting of the Society of Petroleum Engineers of AIME, October 6–9, 1974, Houston, Texas, USA. Society of Petroleum Engineers. [Google Scholar]
 Leverett M. (1941) Capillary behavior in porous solids, Trans. AIME 142, 01, 152–169. [Google Scholar]
 Nguyen V.H., Sheppard A.P., Knackstedt M.A., Pinczewski W.V. (2006) The effect of displacement rate on imbibition relative permeability and residual saturation, J. Pet. Sci. Eng. 52, 1–4, 54–70. [Google Scholar]
 Odeh A.S., Dotson B.J. (1985) A method for reducing the rate effect on oil and water relative permeabilities calculated from dynamic displacement data, J. Pet. Technol. 37, 11, 2–051. [Google Scholar]
 Osoba J.S., Richardson J.G., Kerver J.K., Hafford J.A., Blair P.M. (1951) Laboratory measurements of relative permeability, J. Pet. Technol. 3, 02, 47–56. [Google Scholar]
 Rapoport L.A., Leas W.J. (1953) Properties of linear waterfloods, J. Pet. Technol. 5, 05, 139–148. [Google Scholar]
 Reed J., Maas J. (2018 August) Review of the intercept method for relative permeability correction using a variety of case study data, in: The International Symposium of the Society of Core Analysts, Trondheim, Norway, 26–31 August, 2018. [Google Scholar]
 Richardson J.G., Kerver J.K., Hafford J.A., Osoba J.S. (1952) Laboratory determination of relative permeability, J. Pet. Technol. 4, 08, 187–196. [CrossRef] [Google Scholar]
 Santos I.C.A.B.A., Eler F.M., Nunes D.S.S., Couto P. (2021) Evaluation of capillary end effect in wateroil permeability tests using multiple flow rates technique, Braz. J. Pet. Gas 14, 4, 259–267. [Google Scholar]
 Sidiq H., Amin R., Kennaird T. (2017) The study of relative permeability and residual gas saturation at high pressures and high temperatures, Adv. GeoEnergy Res. 1, 1, 64–68. [Google Scholar]
 Valavanides M.S. (2018) Review of steadystate twophase flow in porous media: independent variables, universal energy efficiency map, critical flow conditions, effective characterization of flow and pore network, Trans. Porous Media 123, 45–99. [Google Scholar]
 Virnovsky G.A., Guo Y., Skjaeveland S.M. (1995 May) Relative permeability and capillary pressure concurrently determined from steadystate flow experiments, in: IOR 19958th European Symposium on Improved Oil Recovery, Vienna, Austria, 15–17 May 1995. European Association of Geoscientists & Engineers, cp107 p. [Google Scholar]
 Virnovsky G.A., Vatne K.O., Skjaeveland S.M., Lohne A. (1998, January) Implementation of multirate technique to measure relative permeabilities accounting, in: SPE Annual Technical Conference and Exhibition, September 27–30, 1998, New Orleans, Louisiana, USA Society of Petroleum Engineers. [Google Scholar]
 Zou S., Liu Y., Cai J., Armstrong R.T. (2020) Influence of capillarity on relative permeability in fractional flows, Water Resour. Res. 56, 11, 1–10. [Google Scholar]
All Tables
Reference input parameters in the simulations. The parameters are based on Kleppe and Morse (1974, January). To scale the capillary pressure, interfacial tension was assumed to be 21 mN/m. The core length was assumed to be 10 cm.
Line analysis of data from Virnovsky et al. (1998) giving corrected saturations and relative permeability points from the intercepts, and slope values C_{s}, C_{o}, C_{w} to be used for derivation of capillary pressure for each fraction F.
System parameters from Virnovsky et al. (1998); Relative permeability correlation parameters to fit the corrected relative permeability data; Jfunction correlation parameters to fit the slope parameters from the same dataset. Experiments marked A, B, C with different velocities reported different permeability.
Line analysis of data from Henderson et al. (1998) giving corrected saturations and relative permeability points from the intercepts, and slope values C_{s}, C_{c}, C_{g} to be used for derivation of capillary pressure for each fraction F. Condensate was less wetting than gas for all fractions, hence C_{g} was not reported for any fractions.
System parameters from Henderson et al. (1998); Relative permeability correlation parameters to fit the corrected relative permeability data; Jfunction correlation parameters to fit the slope parameters from the same dataset.
All Figures
Fig. 1 Illustration of the system including flow and boundary conditions, a typical end effect region and the relevant saturation interval at steady state. denotes the saturation where capillary pressure is zero, while denotes the saturation where the flow function f_{w} equals the injection flow fraction F. 

In the text 
Fig. 2 Input saturation functions based on Kleppe and Morse (1974, January): oil and water relative permeability (left) and scaled and absolute capillary pressure (right). The experimental points are shown together with the correlations based on (84)–(86). 

In the text 
Fig. 3 Fractional flow functions f_{w} for different oil viscosities (the base value of 2.3 cP and increased values 5 and 25 times higher) against water saturation in lin–lin plot (left) and semilog plot (right). 

In the text 
Fig. 4 Validation of the analytical model by comparing results from the analytical solution and a commercial software (Sendra) at high (a) and low (b) injected fraction (F = 0.5 and 0.05, respectively) at different total injection rates, as indicated in the legend, measured in pore volumes per day. 

In the text 
Fig. 5 Saturation profiles Y(s_{w}) calculated using (24) which consists of a dimensionless number N_{0} and a saturation integral. Variation of parameters in N_{0} (a) that increase N_{0} by the same factor from the base case, but do not change the saturation integral, result in an identical new case corresponding to a compression of the base case profile. Changing viscosity (ratio) or flow fraction (b) affect the saturation integral and can change the saturation range and profile shape in general. Compressing these profiles would not cause them to overlap. 

In the text 
Fig. 6 Normalized saturation profiles from the analytical solution (in blue) for different injected flow fractions and oil viscosities and their corresponding approximated saturation profiles assuming polynomial shape for n = 5. The profiles are plotted against Y over the range 0 to 10Y_{av} where Y_{av} and Y_{cee} = (n + 1)Y_{av} are marked with dashed vertical lines left and right, respectively. 

In the text 
Fig. 7 Plots of average saturation vs. inverse total rate (left) and plots of inverse relative permeabilities of oil and water vs. inverse total rate (right) under different flow fractions and oil viscosities. The results are calculated both using the general solutions (50), (54) and (55) which are valid under all conditions, and using the analytical solution (AnSol) in (50) and (80)–(83). For the case F = 0.01, μ_{o} = 2.3 cP simulation results from the commercial (COM) software Sendra are added for comparison. 

In the text 
Fig. 8 Comparison of the general solutions (full or dash/dotted lines) to the analytical solution (straight dashed line) in terms of in (a) and (c) and in (b) and (d) for changes in flow fraction F in (a) and (b) and changes in oil viscosity in (c) and (d). All solutions fall on the same straight line at low dimensionless number (high rates) which scale the data. 

In the text 
Fig. 9 The difference between general solutions and the analytical solution in terms of plotted against C_{s}/N_{0} in (a) and (c) and plotted against in (b) and (d) for changes in flow fraction F in (a) and (b) and changes in oil viscosity in (c) and (d). The general solution begins to deviate from the analytical solution at a range of dimensionless numbers. 

In the text 
Fig. 10 The difference between general solutions and the analytical solution in terms of in (a) and (c) and in (b) and (d) for changes in flow fraction F in (a) and (b) and changes in oil viscosity in (c) and (d). All results are plotted against Y_{cee} and demonstrate that the analytical solution is valid for Y_{cee} < 1. 

In the text 
Fig. 11 Effective relative permeabilities and saturations calculated for fractions F between 0.005 and 0.995 for injection total rates varying from 0.1 to 1000 PV/d. The calculated number N_{0} is shown, where high N_{0} corresponds to high injection rate. Each curve is based on a fixed N_{0}, and each saturation on a given curve is obtained from one value of flow fraction. The reference (“true”) input relative permeabilities and capillary pressure curves were used and were only approximated by the effective curves at sufficiently high N_{0}. End effects are important at low N_{0} and shift the saturations towards = 0.61 and reduce the effective relative permeabilities. 

In the text 
Fig. 12 Interpretation of experimental data (crosses) from Virnovsky et al. (1998) where average saturation (a) and inverse effective relative permeabilities in (b) and (c) and in (d) and (e) are plotted against . Straight lines are drawn through the data for each flow fraction F to find the intercept (marked with circle points) corresponding to saturations and relative permeability points corrected for end effects. The slopes of the lines provide values of C_{s}, C_{o}, C_{w} at each fraction. 

In the text 
Fig. 13 Scaling of the data from Virnovsky et al. (1998) by plotting against , (b) and against (when oil is less wetting) and against (when water is less wetting). 

In the text 
Fig. 14 Relative permeabilities in (a) with points based on intercepts from the steady state measurements from Virnovsky et al. (1998) and the correlations (full lines) using (85) and (86) that best fit the points. The best fitting capillary pressure correlation is shown in (b) as scaled Jfunction and in mbar based on matching the slope parameters from the steady state measurements. The curve is compared to a primary drainage curve measured with centrifuge on the core used in Virnovsky et al. (1998). Our derived Jfunction is more consistent with the experimental data as reflected in a higher residual water saturation and the presence of negative capillary pressures. 

In the text 
Fig. 15 Optimal match of the slope coefficients C_{s}, C_{w}, C_{o} after tuning the Jfunction. 

In the text 
Fig. 16 Interpretation of experimental data (crosses) from Henderson et al. (1998) where average condensate saturation (a) and inverse effective relative permeabilities in (b) and in (c) are plotted against . Straight lines are drawn through the high rate data for each flow fraction F to find the intercept (marked with circle points) corresponding to saturations and relative permeability points corrected for end effects. The slopes of the lines provide values of C_{s}·C_{o}, C_{w} at each fraction. 

In the text 
Fig. 17 Scaling of the data from Henderson et al. (1998) by plotting against , (b) and against (when gas is less wetting) and against (when condensate is less wetting). 

In the text 
Fig. 18 Relative permeabilities in (a) with points based on intercepts from the steady state measurements from Henderson et al. (1998) and the correlations (full lines) using (85) and (86) that best fit the points. The best fitting capillary pressure correlation is shown in (b) as scaled Jfunction and in mbar based on matching the slope parameters from the steady state measurements. 

In the text 
Fig. 19 Optimal match of the slope coefficients C_{s} (a) and C_{c} (b) after tuning the Jfunction. C_{g} was not relevant. 

In the text 