Relative Permeability: A Critical Parameter in Numerical Simulations of Multiphase Flow in Porous Media

: Effective multiphase ﬂow and transport simulations are a critical tool for screening, selection, and operation of geological CO 2 storage sites. The relative permeability curve assumed for these simulations can introduce a large source of uncertainty. It signiﬁcantly impacts forecasts of all aspects of the reservoir simulation, from CO 2 trapping efﬁciency and phase behavior to volumes of oil, water, and gas produced. Careful consideration must be given to this relationship, so a primary goal of this study is to evaluate the impacts on CO 2 -EOR model forecasts of a wide range of relevant relative permeability curves, from near linear to highly curved. The Farnsworth Unit (FWU) is an active CO 2 -EOR operation in the Texas Panhandle and the location of our study site. The Morrow ‘B’ Sandstone, a clastic formation composed of medium to coarse sands, is the target storage formation. Results indicate that uncertainty in the relative permeability curve can impart a signiﬁcant impact on model predictions. Therefore, selecting an appropriate relative permeability curve for the reservoir of interest is critical for CO 2 -EOR model design. If measured laboratory relative permeability data are not available, it must be considered as a signiﬁcant source of uncertainty.


Introduction
With the prospects of climate change looming and an ever-increasing demand for power generation and heavy industry, reducing greenhouse gas emissions from large point-source emitters, such as coal and natural gas power plants or fertilizer operations, has become paramount. Geologic carbon storage (GCS) is one potential path for emission reductions. Carbon dioxide is captured from the large point-source emitters, compressed into a supercritical state and injected into a suitable storage formation such as a depleted oil and gas reservoir or a deep saline aquifer [1][2][3][4][5]. Mature oil fields undergoing CO 2 -Enhanced Oil Recovery (CO 2 -EOR) are another promising option for GCS that also offer an economic benefit. The incremental recovery for CO 2 -EOR operations can produce an additional 7-23% of the oil in place while simultaneously storing roughly 40% of the CO 2 injected [6].
Multiphase flow and transport simulation that can characterize CO 2 effects in oil and water are an integral part of designing GCS projects for oil and gas fields. These simulations are used in project design, permitting, forecasting oil production and storage capacity, and quantifying possible site risks. Understanding the uncertainty in the simulation model inputs and their impact on performance and predictions is critical for project success. The permeability distribution in a reservoir is probably the biggest source of uncertainty, followed closely by relative permeability. Three-phase relative permeability has significant impacts on fluid flow and storage capacity, yet is poorly understood and often generalized. Laboratory measurements have historically been focused on measuring two-and threephase relative permeability curves for oil and gas (CH 4 ) reservoirs [7][8][9]. From these data, empirical models have been developed which promise broad applicability, including use the whole oil field at Farnsworth. Petrophysical properties of porosity and permeability were populated across the domain by stochastic algorithm constrained by well logs [16]. A subset of the domain centered on the west half of the field was used for this study. This is the current active injection and production area and primary focus of the SWP research.

Relative Permeability
For this study, 17 different relative permeability curves were constructed that represent a possible set of curves which may apply to the reservoir. The goal is to study a wide range of relative permeability curves with the assumption that as the relative permeability curve becomes more linear, the saturation end-points and the relative permeability end-points become larger, representing a transition from low fluid mobility to high fluid mobility. The residual phase saturation, maximum relative permeability, and the curve shape were varied to create a range of curves that bracket the parameter space from near linear to highly curved. As of the time of this study, there was only a single relative permeability curve measured for the Morrow Sandstone [21,22]. This lack of measured data necessitates that a wide range of input parameters be examined to understand the influence relative permeability has on a CO 2 -EOR operation. Figure 1 highlights a representative selection of curves that transition from highly curved to near linear. The remaining relative permeability curves used in the study are shown in Figures 2, A1 and A2.
gies 2021, 14, x FOR PEER REVIEW 3 of 13 data collected during the ongoing characterization effort including SEM, XRD, seismic data sets, well logs, core samples, and thin sections [26][27][28][29]. The full geological model represents the whole oil field at Farnsworth. Petrophysical properties of porosity and permeability were populated across the domain by stochastic algorithm constrained by well logs [16]. A subset of the domain centered on the west half of the field was used for this study. This is the current active injection and production area and primary focus of the SWP research.

Relative Permeability
For this study, 17 different relative permeability curves were constructed that represent a possible set of curves which may apply to the reservoir. The goal is to study a wide range of relative permeability curves with the assumption that as the relative permeability curve becomes more linear, the saturation end-points and the relative permeability endpoints become larger, representing a transition from low fluid mobility to high fluid mobility. The residual phase saturation, maximum relative permeability, and the curve shape were varied to create a range of curves that bracket the parameter space from near linear to highly curved. As of the time of this study, there was only a single relative permeability curve measured for the Morrow Sandstone [21,22]. This lack of measured data necessitates that a wide range of input parameters be examined to understand the influence relative permeability has on a CO2-EOR operation. Figure 1 highlights a representative selection of curves that transition from highly curved to near linear. The remaining relative permeability curves used in the study are shown in Figure 2 and Figures A1 and A2. Figure 1. A representative selection of relative permeability curves used in this study. All curves used follow a similar trend of high residual saturation and highly curved to low residual saturation and near linear curve. Only those shown here vary the maximum relative permeability in addition to the other two variables. Figure 1. A representative selection of relative permeability curves used in this study. All curves used follow a similar trend of high residual saturation and highly curved to low residual saturation and near linear curve. Only those shown here vary the maximum relative permeability in addition to the other two variables. In this study we used the Eclipse ® numerical simulation software for all the simulations. A key feature that we leveraged is its ability to use lookup tables with linear interpolation between points instead of empirical formulas to calculate the relative permeability. This effectively decouples the empirical formula from the numerical simulator, allowing much greater flexibility in specifying the saturation and relative permeability endpoints as well as the degree of curvature. Any empirical formula could have been chosen to create the curves.
We chose to use a modified version of the Corey's Curve formula to create all the Eclipse ® lookup tables used except for the Base Case, discussed below. This relationship provides great flexibility for developing a wide variety of curves from highly non-linear to linear with the added ability to modify the residual phase saturation and maximum relative permeability endpoints. The Corey's Curve formula allows the creation of a wide range of relative permeability curves based on three sets of inputs, residual saturation, maximum relative permeability, and lambda (the curve parameter). The gas and nonaqueous liquid (oil) curve are described by Equations (1) and (2), where St is the liquid saturation defined as the non-aqueous liquid (oil) saturation plus the residual aqueous liquid saturation, Sgr is the residual gas saturation, Str is the residual liquid saturation defined as the residual oil saturation plus the residual water saturation, λ defines the shape of the curve, krgmax is the maximum gas relative permeability, krnmax is the maximum non-aqueous liquid (oil) relative permeability, and krg(St) and krn(St) are the gas and nonaqueous liquid relative permeability at liquid saturation St [19].
The water and oil curves are described by Equations (3) and (4), where Sl is the aqueous saturation, Slr is the residual aqueous saturation, Snr is the residual non-aqueous (oil) saturation, λ defines the shape of the curve, krlmax is the maximum aqueous liquid relative permeability, krnmax is the maximum non-aqueous liquid relative permeability, and krl(Sl) and krn(Sl) are the aqueous and non-aqueous liquid relative permeability at aqueous liquid saturation Sl [19]. Tables 1 and 2 have the parameters (residual phase saturation, maximum phase relative permeability, and lambda) used to calculate all the relative permeability curves in this study, including the subset shown in Figure 1. krl(Sl) = krlmax((Sl-Slr)/(1-Slr-Snr)) λ (3) Figure 2. Plot shows the Base Case relative permeability curve from a UNOCAL study (solid lines) and the C1 relative permeability curve that represents the most linear curve used in the study.
In this study we used the Eclipse ® numerical simulation software for all the simulations. A key feature that we leveraged is its ability to use lookup tables with linear interpolation between points instead of empirical formulas to calculate the relative permeability. This effectively decouples the empirical formula from the numerical simulator, allowing much greater flexibility in specifying the saturation and relative permeability endpoints as well as the degree of curvature. Any empirical formula could have been chosen to create the curves.
We chose to use a modified version of the Corey's Curve formula to create all the Eclipse ® lookup tables used except for the Base Case, discussed below. This relationship provides great flexibility for developing a wide variety of curves from highly non-linear to linear with the added ability to modify the residual phase saturation and maximum relative permeability endpoints. The Corey's Curve formula allows the creation of a wide range of relative permeability curves based on three sets of inputs, residual saturation, maximum relative permeability, and lambda (the curve parameter). The gas and non-aqueous liquid (oil) curve are described by Equations (1) and (2), where S t is the liquid saturation defined as the non-aqueous liquid (oil) saturation plus the residual aqueous liquid saturation, S gr is the residual gas saturation, S tr is the residual liquid saturation defined as the residual oil saturation plus the residual water saturation, λ defines the shape of the curve, k rg max is the maximum gas relative permeability, k rn max is the maximum non-aqueous liquid (oil) relative permeability, and k rg (S t ) and k rn (S t ) are the gas and non-aqueous liquid relative permeability at liquid saturation S t [19].
The water and oil curves are described by Equations (3) and (4), where S l is the aqueous saturation, S lr is the residual aqueous saturation, S nr is the residual non-aqueous (oil) saturation, λ defines the shape of the curve, k rl max is the maximum aqueous liquid relative permeability, k rn max is the maximum non-aqueous liquid relative permeability, and k rl (S l ) and k rn (S l ) are the aqueous and non-aqueous liquid relative permeability at aqueous liquid saturation S l [19]. Tables 1 and 2 have the parameters (residual phase saturation, maximum phase relative permeability, and lambda) used to calculate all the relative permeability curves in this study, including the subset shown in Figure 1.
Energies 2021, 14, 2370 The Base Case model permutation uses a relative permeability curve (BC in Figure 2) developed by UNOCAL for a simulation study focused on the efficacy of water and CO 2 flooding at the Farnsworth Unit [21,22]. The UNOCAL study did not provide a capillary pressure curve for the Morrow 'B' Sandstone. For this site, it is initially believed that capillary pressure had a negligible effect on phase movement and was thus ignored in the previous models to aid in simplicity and computational speed. For consistency with previous work and to reduce extra variables that may influence the simulation results capillary pressure was not included [16,25]. We plan to study the influence that the addition of capillary pressure may have on model forecasts to determine if this assumption is valid.

Farnsworth Units Model Domain
The Farnsworth Unit is divided into an east and west half that appear to be hydraulically split [22,23]. The west half is the site of most production and injection operations historically and are where current and future CO 2 -EOR operations are occurring [22,23]. A detailed geological model encompassing the Morrow 'B' Sandstone in the west half of the reservoir was the basis of our simulation model [26]. We up-scaled this geologic model with over 26 million cells to a simulation model consisting of 33,756 active cells. Along with the grid geometry, we up-scaled the permeability and the porosity ( Figure 3). This yielded a mean permeability of 39 mD with a standard deviation of 54 mD and a mean porosity of 14% with a standard deviation of 3%, in line with the original geological model values of 13.6% porosity and 27 mD permeability. We assumed the sealing formation, the Morrow Shale, made a no-flow boundary on all sides as well as the top and bottom of the reservoir, so only the reservoir interval is included in the simulation model. See Moodie et al. (2019) for more details on the model domain description [30]. vious work and to reduce extra variables that may influence the simulation results capillary pressure was not included [16,25]. We plan to study the influence that the addition of capillary pressure may have on model forecasts to determine if this assumption is valid.

Farnsworth Units Model Domain
The Farnsworth Unit is divided into an east and west half that appear to be hydraulically split [22,23]. The west half is the site of most production and injection operations historically and are where current and future CO2-EOR operations are occurring [22,23]. A detailed geological model encompassing the Morrow 'B' Sandstone in the west half of the reservoir was the basis of our simulation model [26]. We up-scaled this geologic model with over 26 million cells to a simulation model consisting of 33,756 active cells. Along with the grid geometry, we up-scaled the permeability and the porosity (Figure 3). This yielded a mean permeability of 39 mD with a standard deviation of 54 mD and a mean porosity of 14% with a standard deviation of 3%, in line with the original geological model values of 13.6% porosity and 27 mD permeability. We assumed the sealing formation, the Morrow Shale, made a no-flow boundary on all sides as well as the top and bottom of the reservoir, so only the reservoir interval is included in the simulation model. See Moodie et al. (2019) for more details on the model domain description [30].
The initial conditions were derived from the results of the history-matched primary production and water-flooding modeling study done previously by the SWP and include the oil saturation and oil component distribution, the water saturation and the pressure [16]. The simulation is initialized with no gas phase. All methane (CH4) and other light volatiles are dissolved in the oil. A compositional fluid model is specified for this study based on a fluid properties report from initial exploration, then refined by fluid modeling [31,32].

Well Operations Schedule and Model Fit
The Farnsworth Unit operates under a water alternating gas (WAG) injection scheme. The schedule used in this study mimics current practices and future plans. Figure 3 indicates the location of the injection and production wells used in the model. The production and injection schedule are broken into two main operation periods. The first period extends from 1 December 2010 to 31 July 2016 and uses historical monthly injection and The initial conditions were derived from the results of the history-matched primary production and water-flooding modeling study done previously by the SWP and include the oil saturation and oil component distribution, the water saturation and the pressure [16]. The simulation is initialized with no gas phase. All methane (CH 4 ) and other light volatiles are dissolved in the oil. A compositional fluid model is specified for this study based on a fluid properties report from initial exploration, then refined by fluid modeling [31,32].

Well Operations Schedule and Model Fit
The Farnsworth Unit operates under a water alternating gas (WAG) injection scheme. The schedule used in this study mimics current practices and future plans. Figure 3 indicates the location of the injection and production wells used in the model. The production and injection schedule are broken into two main operation periods. The first period extends from 1 December 2010 to 31 July 2016 and uses historical monthly injection and production data to define well rates. The second period extends from 1 August 2016 to 1 January 2036 and models potential future operations through to the end of the field's lifetime. During this second period, the CO 2 injection volume gradually decreases until it is completely reliant on recycled CO 2 to meet targets. See

Discussion
Results indicate that the shape and endpoints of the relative permeability curve impart a significant influence on model forecasts. Generally, the more non-linear the curve, the more CO2 that is predicted to be stored, and the less oil, water, and gas is predicted to be produced. However, the saturation endpoints and the assumed maximum relative permeability have a significant impact on this trend, with some of the more non-linear curves (C11 & C12) predicting more oil production than the most linear curve (C1).

Carbon Dioxide Storage
The total amount of stored CO2 forecasted in this modeling study is between 2.8 million tons and 3.2 million tons by the end of the simulation. The largest fraction of the CO2 is dissolved in the oil phase, between 1.2 and 2.4 million tons. The cases using highly nonlinear relative permeability curves with narrow saturation ranges (C11 through C16) had overall a larger portion of the CO2 dissolved in the oil phase while the more linear curves with broader saturation ranges (C1 through C7) have almost a million tons less CO2 dissolved in the oil phase but significantly more CO2 in the supercritical phase. This indicates

Discussion
Results indicate that the shape and endpoints of the relative permeability curve impart a significant influence on model forecasts. Generally, the more non-linear the curve, the more CO 2 that is predicted to be stored, and the less oil, water, and gas is predicted to be produced. However, the saturation endpoints and the assumed maximum relative permeability have a significant impact on this trend, with some of the more non-linear curves (C11 & C12) predicting more oil production than the most linear curve (C1).

Carbon Dioxide Storage
The total amount of stored CO 2 forecasted in this modeling study is between 2.8 million tons and 3.2 million tons by the end of the simulation. The largest fraction of the CO 2 is dissolved in the oil phase, between 1.2 and 2.4 million tons. The cases using highly non-linear relative permeability curves with narrow saturation ranges (C11 through C16) had overall a larger portion of the CO 2 dissolved in the oil phase while the more linear curves with broader saturation ranges (C1 through C7) have almost a million tons less CO 2 dissolved in the oil phase but significantly more CO 2 in the supercritical phase. This indicates that the CO 2 in the supercritical phase has an inverse relationship to the CO 2 dissolved in the oil phase with respect to the changes in the relative permeability curve ( Figure 5). As the curve becomes more non-linear, the fluid mobility decreases, leading to a decrease in the volume of CO 2 dissolved in the oil and a corresponding increase in the supercritical CO 2 . The highest fluid mobility curve (model C1) predicts an almost even CO 2 distribution between the supercritical phase (46%) and oil phase (44%); while the lowest fluid mobility curves (C14 to C16) predict a significant difference in CO 2 phase distribution, 10% in the supercritical phase and 82% dissolved in the oil phase for model C14, see Table 3 for stored CO 2 mass distribution. The CO 2 mass dissolved in the water phase is mostly unaffected by changes in the relative permeability curve, varying by only 2% across all model permutations, from 8% to 10% of the total CO 2 stored. The Base Case (BC) model predicts the largest mass of CO 2 stored in the reservoir and the phase distribution matches the models with relative permeability curves that describe the curve becomes more non-linear, the fluid mobility decreases, leading to a decrease in the volume of CO2 dissolved in the oil and a corresponding increase in the supercritical CO2. The highest fluid mobility curve (model C1) predicts an almost even CO2 distribution between the supercritical phase (46%) and oil phase (44%); while the lowest fluid mobility curves (C14 to C16) predict a significant difference in CO2 phase distribution, 10% in the supercritical phase and 82% dissolved in the oil phase for model C14, see Table 3   . This chart shows CO2 storage for each of the model permutations. 'CO2 SC' is the supercritical CO2, 'CO2 Oil' is the CO2 and CH4 that is dissolved in the oil phase, 'CO2 Water' is the CO2 dissolved in the aqueous phase, and 'CO2 Total' is the total amount of gas, both CO2 and CH4, in the reservoir.  Figure 5. This chart shows CO 2 storage for each of the model permutations. 'CO 2 SC' is the supercritical CO 2 , 'CO 2 Oil' is the CO 2 and CH 4 that is dissolved in the oil phase, 'CO 2 Water' is the CO 2 dissolved in the aqueous phase, and 'CO 2 Total' is the total amount of gas, both CO 2 and CH 4 , in the reservoir.

Oil Production
Oil production does not follow the same trend as the CO 2 storage. Models C1 through C9 indicate a declining oil production as fluid mobility described by the relative permeability curves decreases, except models C4 and C7 ( Table 4). The relative permeability curve used in models C11, C12, and C14 describes a lower fluid mobility condition but predicts similar oil production as the highest fluid mobility models (C1 and C2); whereas the lowest fluid mobility relative permeability models predict the lowest oil production, as would be expected. A possible reason for the high oil production show in models C11, C12, and C14 is that the oil saturation falls within a range on the relative permeability curve that promotes the oil phase mobility, the mid-range of the relative permeability curve. Oil saturations are between 40% and 60% in the active production and injection areas giving relative permeability ranges of 0.23 to 0.35. Within the same area, the relative permeabilities to water never get above 0.1, promoting oil mobility over water mobility. The Base Case (BC) model's oil production falls within the middle of the range predicted by this study, similar to models C6, C10, C13 (Table 4). This indicates that the Base Case relative permeability curve was measured from an area of the field that exhibits medium fluid mobility when compared to the range of fluid mobilities predicted by the synthetic relative permeability curves.

Pressure
The influence of relative permeability on the pressure field was highly time dependent. Figure 6 indicates that during the first phase of injection when historical data are used to control the injection rates (2010 to 2017), there is less difference in pressure between the relative permeability curves tested. During the predictive phase of the injection schedule (2017 to 2036), the difference across all the relative permeability curves tested increased to 27% by the end of injection. This variation in pressure increases significantly as the proportion of recycled CO 2 in the injection stream increases. An inflection point in all of the models on 1 January 2024 marks the point when CO 2 availability for injection becomes tied to production volumes and hence the fluid movement within the reservoir. On 1 January 2030, there is another inflection point that marks when new CO 2 to the model is stopped and only recycled CO 2 is available to the injection wells.
volume of CO2 available to meet injection targets, the impact of the relative permeability curves becomes much more pronounced. Curves that restrict the fluid movement through a high degree of curvature and narrow saturation range, such as C15 and C16, have lower oil production, less CO2 present in the reservoir, and a higher average pressure, while the curves that promote fluid movement through more linear curvature and a wider saturation range (C1 and C2) have higher oil production, more CO2 in the reservoir, and a lower average reservoir pressure.

Conclusions
Results of this study indicate that small variations in the shape of the relative permeability curve have a significant impact on the model forecasts. While all models predicted nearly the same total CO2 stored in the reservoir, the phase it is stored as (supercritical vs. dissolved in oil vs. dissolved in water) is greatly influenced by the relative permeability curve. Relative permeability curves that describe low fluid mobility predict most of the CO2 is stored in the oil phase with very little in the supercritical phase. The relative permeability curves that describe the highest overall fluid mobility predicts that there is an even distribution of CO2 in the supercritical phase and the oil phase, allowing the CO2 to migrate faster and in greater quantities to the production wells, leading to lower amounts of total stored CO2. The higher the mobility, the more contact between the CO2 plume and the oil, increasing the amount of CO2 that is dissolved in the oil phase and thereby increasing production when compared to the relative permeability curves that predict low overall fluid mobility. The reduction in the oil's viscosity due to CO2 dissolution allows it The relative permeability curves appear to have a smaller, but still significant, impact on the average reservoir pressure when there is an unlimited source of CO 2 for the injection wells to meet their targets. Once the fluid mobility within the reservoir influences the volume of CO 2 available to meet injection targets, the impact of the relative permeability curves becomes much more pronounced. Curves that restrict the fluid movement through a high degree of curvature and narrow saturation range, such as C15 and C16, have lower oil production, less CO 2 present in the reservoir, and a higher average pressure, while the curves that promote fluid movement through more linear curvature and a wider saturation range (C1 and C2) have higher oil production, more CO 2 in the reservoir, and a lower average reservoir pressure.

Conclusions
Results of this study indicate that small variations in the shape of the relative permeability curve have a significant impact on the model forecasts. While all models predicted nearly the same total CO 2 stored in the reservoir, the phase it is stored as (supercritical vs. dissolved in oil vs. dissolved in water) is greatly influenced by the relative permeability curve. Relative permeability curves that describe low fluid mobility predict most of the CO 2 is stored in the oil phase with very little in the supercritical phase. The relative permeability curves that describe the highest overall fluid mobility predicts that there is an even distribution of CO 2 in the supercritical phase and the oil phase, allowing the CO 2 to migrate faster and in greater quantities to the production wells, leading to lower amounts of total stored CO 2 . The higher the mobility, the more contact between the CO 2 plume and the oil, increasing the amount of CO 2 that is dissolved in the oil phase and thereby increasing production when compared to the relative permeability curves that predict low overall fluid mobility. The reduction in the oil's viscosity due to CO 2 dissolution allows it greater fluid mobility and may be why there is an increase in oil production when compared to some of the relative permeability curves that predict low fluid mobility. This may also account for why there is less CO 2 stored in the oil phase with the relative permeability curves that describe high fluid mobility.
The findings of this study indicate that the relative permeability curve is a critical parameter that must be given careful consideration when designing multiphase flow models. It is essential to understand and quantify the uncertainty in the relative permeability curve.
If measured laboratory relative permeability data are not available or limited for the study domain, the relative permeability curve should be considered a significant source of model uncertainty and accounted for as part of the simulation effort.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Appendix A
greater fluid mobility and may be why there is an increase in oil production when compared to some of the relative permeability curves that predict low fluid mobility. This may also account for why there is less CO2 stored in the oil phase with the relative permeability curves that describe high fluid mobility.
The findings of this study indicate that the relative permeability curve is a critical parameter that must be given careful consideration when designing multiphase flow models. It is essential to understand and quantify the uncertainty in the relative permeability curve. If measured laboratory relative permeability data are not available or limited for the study domain, the relative permeability curve should be considered a significant source of model uncertainty and accounted for as part of the simulation effort.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.
Appendix A Figure A1. The C3, C6, C9, C12, and C15 relative permeability curves. The saturation endpoints and curvature are varied, while the relative permeability end points remain constant. Figure A1. The C3, C6, C9, C12, and C15 relative permeability curves. The saturation endpoints and curvature are varied, while the relative permeability end points remain constant.