Mathematical modeling of hot carbonated waterflooding as an enhanced oil recovery technique

International Journal of Multiphase Flow 115 (2019) 181–195

Contents lists available at ScienceDirect

International Journal of Multiphase Flow journal homepage: www.elsevier.com/locate/ijmulflow

Mathematical modeling of hot carbonated waterflooding as an enhanced oil recovery technique Tamires A. Marotto, Adolfo P. Pires∗ Laboratório de Engenharia e Exploração de Petróleo - Universidade Estadual do Norte Fluminense Darcy Ribeiro, Macaé, RJ, Brazil

a r t i c l e

i n f o

Article history: Received 15 October 2018 Revised 22 March 2019 Accepted 29 March 2019 Available online 3 April 2019 Keywords: Enhanced oil recovery Thermal methods Miscible flooding Hot carbonated waterflooding Solvent and hot water co-injection

a b s t r a c t Due to environmental concerns, carbon dioxide has been increasingly used in Enhanced Oil Recovery (EOR) projects. Among other options, it can be dissolved in injected water and due to phase equilibrium conditions, part of the injected carbon dioxide is exchanged with the oil phase, creating a miscible EOR process. Thermal EOR methods are the most suitable, and sometimes the only option to produce viscous reservoir fluids. Recently, combined EOR methods have called the attention of oil companies and researchers. In this work we present the analytical solution for the problem of oil displacement by hot carbonated waterflooding, a combined thermal-miscible EOR technique. We consider one-dimensional, two-phase, three-component (oil, carbon dioxide, and water) flow in a homogeneous and isotropic porous medium. Other hypotheses of the model are incompressible system with no diffusion and no chemical reactions; gravitational, dispersive and capillary effects are neglected. Following Amagat’s law, we do not consider volume of mixing, and Henry’s law is used to model the solvent (carbon dioxide) distribution between phases. The heat capacities of the components and the rock were considered constant. The dependent variables of the problem are oil saturation, carbon dioxide concentration in the oil phase, and temperature. The solution of this 3 × 3 quasi-linear hyperbolic system is composed of shock and rarefaction waves and constant states, and it was obtained using the method of characteristics. Solutions for different relations between Henry’s constants were developed, and a sensitivity analysis for the CO2 concentration was performed. It is shown that there is an optimum parameter (relation between Henry’s constants) for obtaining the highest recovery factor. The efficiency of this technique is compared to hot waterflooding. © 2019 Elsevier Ltd. All rights reserved.

1. Introduction World oil reserves are decreasing and discoveries of new fields are becoming rare and costly due to the geologic complexity of new frontiers (Zhao, 2004). Under this scenario, it is fundamental to increase the recovery factor of producing fields (Sohrabi et al., 2012). Enhanced Oil Recovery (EOR) processes have been applied in reservoirs since the 1940’s, due to the significant amounts of remaining oil reserves after primary and secondary recovery (Green and Willhite, 1998). Since then, interest in EOR has grown owing to the increase in oil reserves and production. Enhanced Oil Recovery can be defined as a group of techniques used to improve the recovery of hydrocarbons by the injection of materials that are not normally present in the reservoir (Lake, 1989). Enhanced Oil Recovery may also be understood as a technique that modifies porous medium interactions with the flowing phases, changing the

∗

Corresponding author. E-mail address: [email protected] (A.P. Pires).

https://doi.org/10.1016/j.ijmultiphaseflow.2019.03.024 0301-9322/© 2019 Elsevier Ltd. All rights reserved.

relative permeabilities or phase viscosities, thereby raising the oil recovery. Most EOR methods may be classified into three main groups; chemical, thermal and miscible methods. Chemical flooding improves oil recovery by lowering the interfacial tension between oil and water, or reducing the water mobility (Ezekwe, 2011). Injection of a hot fluid (thermal methods) decreases the viscosity of reservoir oil through heat exchange. For several reservoirs, it is the only way to produce its viscous oil. Thermal methods of EOR are responsible for more than half of the total EOR oil production (Kokal and Al-Kaabi, 2010). These methods include steam injection, hot water injection and in-situ combustion. Most thermal oil production comes from cyclic steam injection and steamdrive, which are defined as steam injection methods (Green and Willhite, 1998). Miscible methods are based on the injection of a solvent that changes capillary and interfacial forces through mass transfer between displacing and displaced phases. The efficiency of this process is highly dependent on the phase equilibrium. Miscible methods include carbon dioxide injection, continuously or as

182

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

Nomenclature List of symbols A− Upstream shock state A+ Downstream shock state cij Mass fraction of component i in phase j cˆi j Volume fraction of component i in phase j Cpi Heat capacity of component i Cpr Heat capacity of the rock D Shock velocity E Parameter used in Eq. (37) F Parameter used in Eqs. (51) and (56) fj Fractional flow of phase j G Parameter used in Eq. (56) H Parameter used in Eqs. (51) and (56) Hj Enthalpy of phase j Hr Rock enthalpy I Parameter used in Eqs. (51) and (56) J Parameter used in Eq. (37) k Porous medium permeability krj Relative permeability of phase j k0r j Endpoint relative permeability of phase j KD Ratio between the values of Henry’ s constants Ks,j Henry’s constant for solvent s in phase j L Reservoir length nc Number of components nj Corey’s exponent np Number of phases P System pressure r Right eigenvector sj Saturation of phase j s∗j Normalized saturation of phase j srj Residual saturation of phase j T Temperature t Time coordinate tD Dimensionless time variable uj Velocity of phase j uT Total velocity x Space coordinate xD Dimensionless spatial variable Greek characters α (i) Multiplier λ Characteristic velocities (eigenvalues) μj Viscosity of phase j ξ Dimensionless velocity ρj Density of phase j ρˆi Pure component density at P and T ρr Rock density φ Porosity

the gas part of water-alternated-gas injection process. This technique has become more important recently due to carbon dioxide sequestration. Mathematically, EOR methods are modeled by systems of partial differential equations composed by dissipative and convective terms. When dissipative effects due to capillarity, compressibility, and molecular diffusion are neglected, only convective transport is taken into account. So, the mathematical problem is described by a system of non-linear, hyperbolic, first-order differential equations. If the injection and initial conditions are constant and uniform, it is a Riemann problem, and one of the most popular techniques to solve these equations is the method of characteristics (MOC) (Dindoruk, 1992; Seto, 2007).

The classical waterflooding problem, where oil is displaced by continuous injection of water into the reservoir, considers onedimensional and isothermal flow, with no capillary and gravity terms (Buckley and Leverett, 1942). To find the correct saturation profile for an S-shape fractional flow function, a shock (discontinuity) obeying Rankine–Hugoniot conditions is included in the solution obtained for uniform initial water saturation. Similar problems are known as Buckley–Leverett type in reservoir engineering literature. Based on the Buckley–Leverett solution, Welge (1952) developed a graphical method for computing the average water saturation, and hence the oil recovery. The method of characteristics was used to model surfactant flooding in oil reservoirs (Fayers and Perrine, 1958). The mathematical system was based on Darcy’s law and on mass conservation equations. For this chemical flooding model it was noted that besides a Buckley–Leverett type saturation shock, a concentration shock also occurs. So, there is a discontinuity in the saturation profile coinciding with a discontinuity in detergent concentration profile. In this problem, adsorption was considered a function of concentration, furthermore, a concentration dependent relative permeability function was incorporated into the model. The solution of the energy balance equation is the mathematical tool to estimate the thermal behavior of hot fluids injection processes (Spillette, 1965). The mechanisms of heat transfer that may be included are energy transfer due to the physical movement of the injected fluids; thermal conduction from the warmer to the colder regions of the system; and convective transfer between the injected fluid and the original reservoir fluid and rock. The reservoir response to hot water injection described by a two-dimensional mathematical model, composed by mass-balance equations and energy-balance equation, showed that hot waterflooding could improve the displacement efficiency if compared to cold waterflooding (Spillette and Nielsen, 1968). The graphical method developed by Welge was extended to calculate miscible displacement of oil by gas or carbonated water, considering the effect of mass transfer between phases. The carbonated waterflooding case involves mass transfer of CO2 to the oil phase, which changes oil-phase mobility and fractional flow curve (Claridge and Bondor, 1974). Using the fractional flow theory, the problem of oil displacement by carbonated waterflooding in a homogenous and isothermal porous medium was solved by the method of characteristics. Based on mass conservation laws it was possible to calculate the front velocity of water saturated with CO2 and the aqueous phase saturation velocity (Pope, 1980). Many researchers also developed solutions for displacement problems in porous media considering variable temperature. In other words, the reservoir temperature may not be the same as that of the injected fluid. For models involving temperature changes, the phase behavior is affected because the thermodynamic equilibrium depends on temperature. Temperature also affects multiphase flow by changing physical properties of the phases, such as densities and viscosities (Zhu, 2003). Hot waterflooding has been used for heavy oil recovery since the 60’s due to viscosity reduction and thermal swelling of the reservoir oil, and even decreasing the residual oil saturation (Willman et al., 1961; Martin et al., 1968; Holke and Huebner, 1971). The displacement of oil by high flow rate hot water injection in a homogeneous isotropic porous medium was modeled by a system composed of two first-order non-linear equations. The solution of this problem was found by the method of characteristics for cases where thermal capacity is constant or depends on the temperature. Besides the Buckley–Leverett type shock, a temperature shock was also observed (Fayers, 1962). Applying the fractional flow theory, Pope (1980) solved the problem of oil displacement by hot waterflooding in a homogenous porous medium by the method of characteristics. The velocity

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

of the temperature wave was obtained coupling the mass conservation and energy conservation equations. The same problem was solved graphically based on a different set of variables, which depends on specific heat, density and porosity (Lake, 1989). The system of conservation laws governing steam and nitrogen injection into a one-dimensional porous medium saturated by water was proposed by Bruining and Marchesin (2006). It was based on the following hypothesis: incompressible system, no heat transfer by conduction and no capillary effects. In this work, the condensation front and shock waves were analyzed, and it was found that there are four possible types of solutions for a specific set of initial and boundary conditions. Further, Lambert et al. (2010) evaluated all possible solution structures for any initial and boundary conditions. A method to find fundamental solutions for thermal compositional flows in porous media was also developed. In addition, an expression for the evaporation rarefaction wave was derived. This kind of wave occurs in two-phase regions containing liquid water and gas (mixture of nitrogen and steam). A more general model for the nonisothermal Buckley–Leverett problem including a tracer, that may be used to track the flood fronts, was developed for cold and hot waterflooding. The solution presents two shocks: a Buckley–Leverett type and a tracer concentration shock, which is a function of the number of injected porous volumes and the initial water saturation in the transition zone. This problem was solved analytically using the method of characteristics, for both radial and linear geometry (Dindoruk and Dindoruk, 2008). Different combined thermal and miscible methods, such as the co-injection of solvent with hot water, have been proposed trying to take advantage of the thermal (hot waterflooding) and solvent methods of EOR. Hot carbonated waterflooding becomes a natural option for this type of combined EOR technique. Field applications of carbonated water injection (CWI) have been reported since the 60’s in Oklahoma and Texas, and the results showed that oil recovery was around 50% higher when compared to waterflooding (Christensen, 1961; Hickok and Ramsey Jr., 1962). Injection of carbonated water and water to evaluate the imbibition mechanisms in limestone core samples confirmed that dissolution of CO2 in water improves both oil production rate and ultimate oil recovery. The results were found to be a function of the amount of CO2 dissolved in the imbibition water (Perez et al., 1992). The HOT − CO2 process combines thermal and solvent EOR techniques. In this method, CO2 is heated above the reservoir temperature and injected in the reservoir to reduce oil viscosity and improve its mobility. Application of this technique into a heavy oil reservoir in Bati Raman field, Turkey, decreased oil viscosity and the dissolution of CO2 in reservoir oil led to a higher recovery factor when compared to CO2 flooding (Picha, 2007). An integrated experimental and theoretical study on the application of CWI showed higher oil recovery compared to water injection. The additional oil recovered was obtained due to oil swelling, viscosity reduction, and coalescence of the isolated oil as a result of CO2 diffusion (Riazi et al., 2009). High-pressure laboratory flow experiments found that CWI, if compared to conventional water injection, improves oil recovery both as a secondary or as a tertiary (after waterflooding) recovery method. The main mechanisms contributing to the additional oil production were oil swelling and viscosity reduction due to CO2 diffusion from the carbonated water into the oil phase (Riazi et al., 2011a; 2011b). Coreflood experiments in sandstone samples showed that CWI yields higher oil recovery than conventional waterflooding. However, secondary CWI led to greater oil recovery than tertiary CWI. Besides this, CWI may be used as a CO2 storage injection strategy (Kechut et al., 2011). Similar experiments in cores saturated with

183

light oil achieved a recovery factor of 33%, while the recovery factor found in cores saturated with viscous oil was 12%. This result is mainly attributed to the light oil swelling and coalescence when compared to the viscous oil behavior (Sohrabi et al., 2011). Carbonated waterflooding experiments using three different cores and three different reservoir fluids (a light oil (n-decane), a refined viscous oil and a stock-tank crude oil) were carried out to evaluate the performance of CWI and to quantify the additional oil recovery and CO2 storage. The ultimate oil recovery obtained by CWI, as secondary and tertiary method, is higher than the recovery factor obtained by conventional waterflooding. The displacement of light oil by CWI was also more efficient than the displacement of viscous oil (Sohrabi et al., 2012). Mosavat and Torabi (2014) carried out light oil secondary carbonated waterflooding in sand-packs. Compared to conventional waterflooding, CWI led to a higher oil recovery mainly due to CO2 mass transfer from carbonated water to the oil phase. At the same pressure conditions, the recovery factor for CWI was 79% of Original Oil in Place (OOIP) and for conventional waterflooding it was 60%. The amount of oil recovered and CO2 stored during CWI depend on temperature. Carbonated water injection laboratory experiments carried out in an Iranian carbonate reservoir core showed that the recovery factor using this technique as a secondary or as a tertiary recovery process was higher than the recovery factor of waterflooding. Besides that, the CO2 storage capacity for both tertiary and secondary CWI were almost the same (Shakiba et al., 2016). In this work we propose the use of carbon dioxide as the solvent in the co-injection of solvent and hot water as a combined miscible-thermal EOR process. Dissolved carbon dioxide swells reservoir oil and decreases its viscosity and density. It is also highly soluble in water and may increase formation permeability (Nasehi and Asghari, 2010). In addition, carbon dioxide as an EOR solvent is environmental friendly because it reduces green house gas emission (Dong et al., 2011). The mathematical problem is a system of three hyperbolic equations representing the conservation of oil phase volume, solvent mass and energy. The hyperbolic system was solved using the method of characteristics (Polyanin and Zaitsev, 2003), and different structures of the solution are presented. The next section describes the mathematical model, followed by the solution for a particular set of reservoir-fluid properties. Then, a sensitivity analysis is performed, addressing the effect of the amount of dissolved CO2 in the injected water, and the recovery factor of the proposed technique is compared to hot waterflooding. Finally, conclusion are presented. 2. Physical and mathematical model description The system of governing equations that models the injection of a hot fluid containing a solvent into an oil reservoir consists of oil, solvent and water mass balance and energy conservation. We assumed the following hypotheses: •

• • • • •

• •

•

One-dimensional two-phase flow in a homogeneous porous medium; No diffusion, no chemical reactions; Incompressible system; Gravity and capillary effects are neglected; Local thermal equilibrium; Enthalpies are functions of phases compositions and temperature; Constant heat capacity; Viscosities are functions of phases compositions and temperature; Pure component density is the same in all phases;

184 • •

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

Amagat’s law is valid; There is mass transfer of solvent component between the phases.

Under these assumptions, the mass conservation of components oil, solvent and water can be written as:

∂ ∂ (φρocooso ) + (ρocoouo ) = 0, ∂t ∂x

(1)

∂ ∂ [φ (ρw csw sw + ρocsoso )] + (ρ c u + ρocsouo ) = 0 ∂t ∂ x w sw w

(2)

and

∂ ∂ (φρw cww sw ) + (ρw cww uw ) = 0 ∂t ∂x

(3)

where φ is the porosity, ρ j is the density of phase j, cij is the mass fraction of component i in phase j, sj is the saturation of phase j, uj is the velocity of phase j, t is the time coordinate and x is the space coordinate. The energy conservation is given by:

∂ [φ (ρosoHo + ρw sw Hw ) + (1 − φ )ρr Hr ] ∂t ∂ + (ρ H u + ρw Hw uw ) = 0, ∂x o o o

(4)

Hw = (C pw cww + C ps csw )T ,

(5)

Ho = (C pocoo + C ps cso )T

(6)

and

(7)

where Cpi is the heat capacity of component i, Cpr is the heat capacity of the rock, and T is the temperature. We define the dimensionless time (tD ) and spatial (xD ) variables as:

t 0

u T ( τ )d τ and xD = xL . φL

(8)

The fractional flow function of phase j is defined as the ratio between the velocity of phase j and total velocity:

uj , uT

(9)

where fj is the fractional flow of phase j and uT is the total velocity. Darcy’s law for one-dimensional multi-phase flow neglecting gravity effects is given by:

u j = −k

kr j ∂ P , μj ∂x

(10)

where k is the porous medium permeability, krj is the relative permeability of phase j, and μj is the viscosity of phase j. Considering that Amagat’s Law (Prausnitz et al., 1986) is valid, and that the pure component density is the same in all phases, we can replace the mass fraction by the volume fraction of component i in phase j in the conservation laws using the following equation:

cˆi j =

(12)

where Ks,j is the Henry’s constant for solvent s in phase j. Henry’s constant of carbon dioxide in both phases was determined using Harvey’s model (Harvey, 1996). For oil phase it was calculated at n-hexadecane saturation pressure at system temperature corrected to system pressure applying Poynting’s correction (Prausnitz et al., 1986), and for water phase it is calculated at water saturation pressure at system temperature corrected to system pressure applying Poynting’s correction (Prausnitz et al., 1986). System pressure is defined as the average pressure between injection pressure and initial reservoir pressure, and system temperature was defined the same way. Water saturation pressure and properties can be found in Saul and Wagner (1987). We also use the following auxiliary relations: np np nc    cˆi j = 1, s j = 1 and f j = 1. j=1

(13)

j=1

Applying Eqs. (5)–(11) in Eqs. (1)–(4), we obtain the following system for oil, solvent and water mass balance and energy conservation:

  ∂  ∂  cˆooso + cˆoo fo = 0, ∂ tD ∂ xD

(14)

  ∂  ∂  cˆ s + cˆ s + cˆ f + cˆ f = 0, ∂ tD sw w so o ∂ xD sw w so o

(15)

  ∂  ∂  cˆ s + cˆ f = 0 ∂ tD ww w ∂ xD ww w

(16)

and

Hr = C pr T ,

fj =

cˆsoKs,o = cˆsw Ks,w ,

i=1

where Hj is the enthalpy of phase j, ρ r is the rock density and Hr is the rock enthalpy. The enthalpies of water, oil and rock are defined as:

tD =

Henry’s Law (Prausnitz et al., 1986) to calculate the fugacity of the solvent in the liquid phases (oil and water), we get:

ci j ρ j , ρˆi

   (1 − φ ) ∂T  MTocˆooso + MTs cˆsoso + MTs cˆsw sw + MTw cˆww sw + MTr ∂ tD φ  ∂T  + M cˆ f + MTs cˆso fo + MTs cˆsw fw + MT w cˆww fw = 0, (17) ∂ xD To oo o

where

MT i = ρˆiC pi .

Using Eqs. (12) and (13) in Eqs. (14)–(16) and expanding, we derive the following system for oil, solvent and water mass balance:

 ∂ fo  ∂ so ∂  ∂  − socˆso + − focˆso = 0, ∂ tD ∂ tD ∂ xD ∂ xD   ∂  cˆ s ( 1 − so ) + ∂ tD so o    ∂ Ks,o ∂  + cˆso cˆ f = 0 (1 − f o ) + ∂ xD Ks,w ∂ xD so o

∂ ∂ tD

and



Ks,o cˆso Ks,w

(19)





(20)

 ( 1 − so ) +   ∂ ∂ Ks,o cˆso (1 − f o ) − (1 − fo ) = 0. (21) ∂ xD ∂ xD Ks,w

∂ ∂ ( 1 − so ) − ∂ tD ∂ tD

(11)

where cˆi j is the volume fraction and ρˆi is the density of pure component i at P and T. To relate the solvent concentrations in water and oil phases we consider infinity dilution model for both phases. Therefore, using

(18)

Ks,o cˆso Ks,w



The oil (Eq. (14)) and solvent (Eq. (15)) equations can be written as:

  ∂ f o ∂ so ∂  ∂  socˆso + focˆso = + ∂ tD ∂ xD ∂ xD ∂ tD

(22)

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

and



   ∂ Ks,o ∂ Ks,o − cˆso cˆso ( 1 − so ) − (1 − f o ) ∂ tD Ks,w ∂ xD Ks,w   ∂  ∂  = cˆ s + cˆ f . ∂ tD so o ∂ xD so o

utD + AuxD = 0,



where



u = so cˆso T (23)

Replacing Eqs. (22) and (23) in Eq. (21) shows that the water conservation law is a linear combination of the oil and solvent equations and does not need to be solved. At this point we emphasize that it is necessary to determine   only one concentration cˆso , because cˆsg can be calculated through Henry’s Law (Eq. (12)) and the other concentrations can be found from the first auxiliary relation (Eq. (13)) for oil and water phases. Thus, the unknowns of this problem are oil saturation, solvent concentration in oil phase and temperature. So, the final system is composed by Eqs. (19), (20) and (17):

  ∂  ∂  s 1 − cˆso + f 1 − cˆ = 0, (24) ∂ tD  o ∂ xD o  so   ∂ Ks,o ∂ Ks,o cˆ cˆ (1 − so ) + cˆsoso + (1 − fo ) + cˆso fo = 0 ∂ tD so Ks,w ∂ xD so Ks,w

T

∂ ∂ tD







so 1 − cˆso T MTo + socˆsoT MTs + (1 − so )cˆso

λ (1 )

  1 ∂ f o (1 ) = , r = 0 , ∂ so 0

Applying the chain rule in Eqs. (24)–(26), we obtain

    ∂ so ∂ cˆso ∂ f o ∂ so Ks,o Ks,o + + 1− cˆso ( 1 − so ) ∂ tD Ks,w ∂ tD Ks,w ∂ so ∂ xD   ∂ fo ∂ cˆso Ks,o Ks,o + 1− cˆso + (1 − f o ) Ks,w ∂ xD ∂ cˆso Ks,w   ∂ fo ∂ T Ks,o + 1− cˆso =0 Ks,w ∂ T ∂ xD

and

Ks,o cˆso Ks,w



⎡ cˆso

⎤

Ks,o ( Ms − Mw ) + Mw Ks,w

⎢ ⎥   ⎥ ⎢ fo +  ⎢ ⎥ Ks,o ⎢ ⎥ 1 + cˆso M − M + M − 1 − M ( w s) s w ⎥ ∂T Ks,w ∂T ⎢ ⎥ +⎢ ⎥ ∂ xD = 0. Ks,o ∂ tD ⎢ (1 − φ ) ⎢ ⎥ cˆso Ms − Mw ) + Mw + Mr ( ⎢ ⎥ Ks,w φ ⎢ so +    ⎥ ⎣ ⎦ Ks,o 1 + cˆso ( Mw − Ms ) + Ms − 1 − Mw Ks,w

Rewriting system (Eq. (27)) in the following form:



⎡

⎤ (30)

and

fo +



 1 + cˆso

so +



cˆso

Ks,o ( Ms − Mw ) + Mw Ks,w



Ks,o ( Mw − Ms ) + Ms − 1 − Mw Ks,w

Ks,o (1 − φ ) Mr ( Ms − Mw ) + Mw + Ks,w φ



1 + cˆso

⎡

cˆso

⎤



Ks,o ( Mw − Ms ) + Ms − 1 − Mw Ks,w

1 0 ⎢ ⎥ r (3 ) = ⎢ λ (3 ) − λ (1 ) ⎥ ⎣ ⎦ ∂ fo ∂T

 ,



(31)

where

Ks,o ( f o − so ) Ks,w C=− . Ks,o  so + 1 − so − cˆso Ks,w

    ∂ so  ∂ f o ∂ so ∂ cˆso 1 − cˆso − so + 1 − cˆso ∂ tD ∂ tD ∂ so ∂ xD        ∂ fo  ∂ fo ∂ T ∂ cˆso + 1 − cˆso − fo + 1 − cˆso = 0, ∂ xD ∂ T ∂ xD ∂ cˆso 1−

(29)

1  Ks,o  fo + 1 − fo − cˆso ⎢ λ (2 ) − λ (1 ) ⎥ Ks,w ⎢ ⎥ (2 ) λ (2 ) = ∂ fo ⎥ , r = ⎢ Ks,o  + C ⎣ ⎦ so + 1 − so − cˆso ∂ cˆso Ks,w 0

Ks,o T MTs Ks,w

Ks,w



,

λ (3 ) =

  Ks,o (1 − φ ) + (1 − so ) 1 − cˆso T MT w + T MT r Ks,w φ    ∂ Ks,o + f 1 − cˆso T MTo + focˆsoT MTs + (1 − fo )cˆso T MTs ∂ xD o Ks,w   Ks,o + (1 − fo ) 1 − cˆso T MT w = 0 . (26)



(28)

we obtain an upper triangular matrix, so the characteristic velocities (eigenvalues) of this system are the elements of the main diagonal. Thus, the corresponding eigenpairs are given by:

(25) and

185

(32)

From the eigenpairs, we obtain the multipliers (Bedrikovetsky, 1993) to calculate the rarefaction waves. So, from the first eigenpair, λ(1) and r(1) we get:



α (1 ) =

∂ 2 fo ∂ s2o

−1

(33)

and the rarefactions of the first family are given by

dso = α (1 ) , dξ

(34)

dcˆso =0 dξ

(35)

and

(27)

dT = 0, dξ

(36)

xD . Thus, cˆso and T are constant along the characteristD tics of the first family. For the second eigenpair, we have: where ξ =

∂ fo J +C ∂ cˆso  ,   α (2 ) =     Ks,o  (1 ) Ks,o ∂ fo Ks,o ∂ fo +C 1− J λ − E + λ (2 ) − λ (1 ) 1− J + (E − J ) Ks,w Ks,w ∂ cˆso Ks,w ∂ cˆso 2

(37)

186

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195



where

 Ks,o  1 − so − cˆso Ks,w

J = so + and

E = fo +

(38)

 Ks,o  1 − fo − cˆso . Ks,w

Ks,o 1− Ks,w



Ks,o 1− Ks,w



 fo+

 −

− cˆso

−

− cˆso

 s+ o



Ks,o + Ks,w Ks,o + Ks,w



Ks,o 1− Ks,w



Ks,o 1− Ks,w



 fo−

 s− o

 λ (2 ) − λ (1 )        s,o + C 1 − KKs,w J λ (1 ) − E + λ (2 ) − λ (1 ) 1 − J2

(50) and

(40)

∂ cˆso

Ks,o + Ks,w



(39)

dso = α (2 ) , dξ

 ∂ fo



D= + cˆso

The second family of rarefaction waves is given by

dcˆso = dξ

Ks,o + Ks,w

+ cˆso

+  −  +  − cˆsoF + H fo+ − cˆso F + H fo− + cˆso − cˆso I       D= , + + − − + −

(51)

cˆsoF + H so − cˆsoF + H so + cˆso − cˆso I



Ks,o Ks,w



J ∂∂cˆfsoo +

Ks,o Ks,w

 (E − J )

(41)

and

dT = 0. dξ

(42)

Therefore, T is constant along the characteristics of the second family. From the last eigenpair, λ(3) and r(3) , the multiplier is

α (3 )

= 0.

(43)

In this case, ∇λ(3 ) · r (3 ) = 0 for all so , cˆso and T, a linearly degenerate wave (LeVeque, 2002). From the above mentioned rarefaction waves, we conclude that oil saturation changes while solvent concentration and temperature remain constant along the first rarefaction wave. For the second rarefaction family, both oil saturation and solvent concentration change, and temperature is constant. The oil saturation and temperature change while the solvent concentration is constant along the third family, a degenerate wave (shock). The Rankine–Hugoniot (RH) conditions describe the relationship between the states of the unknowns on both sides of a shock wave. Theses equations are derived from the conservation laws in integral form (LeVeque, 2002). So, the RH condition for this problem is:









socˆoo D = focˆoo ,





F=



Ks,o (MT w − MTs ) + MTs − MTo , Ks,w



(53)

and

I=

Ks,o (MTs − MT w ). Ks,w

(54)

Finally, for the degenerate (Eqs. (44)–(46)) are given by:

D=

wave,

the

jump

conditions

fo+ − fo− − s+ o − so

(55)

and D=

T+





 

cˆsoF + H fo+ + cˆsoI + MTw − T −

T+





− cˆsoF + H s+ o + cˆso I + G − T



  cˆsoF + H fo− + cˆsoI + MTw 

−  , cˆsoF + H so + cˆsoI + G

(56) where

G = MT w +



(52)

H = MTo − MT w

(44)

sw cˆsw + socˆso D = fw cˆsw + focˆso

(45)

(1 − φ ) MT r . φ

(57)

3. Solution

and





soMTocˆooT + soMTs cˆsoT + sw MTs cˆsw T + sw MT w cˆww T + =



(1 − φ ) MTr T D φ 

foMTocˆooT + foMTs cˆsoT + fw MTs cˆsw T + fw MTw cˆww T ,

(46)

where [A] = A+ − A− represents the jump of A from a downstream condition to an upstream condition. For example, the RH condition for Eq. (44) is



where







+ − − + + − − s+ o cˆoo − so cˆoo D = f o cˆoo − f o cˆoo .

(47)

For the first shock type, solvent concentration and temperature + − are constant, therefore cˆso = cˆso = cˆso and T + = T − = T . So, keeping solvent concentration in oil phase and temperature constant in Eqs. (44)–(46), the jump condition for this case reduces to:

D=

fo+ s+ o

− −

fo− . s− o

(48)

Following the same procedure for the second type shock, where temperature is constant, leads to the following system of equations:

  +   − + − 1 − cˆso fo − 1 − cˆso f    o , D=  + + − − 1 − cˆso so − 1 − cˆso so

(49)

In this section we present different solutions for this problem. Corey’s model (Corey et al., 1956) was used to calculate the relative permeability of the phases:

 n j

kr j = k0r j s∗j

,

(58)

where krj is the relative permeability of phase j, k0r j is the endpoint relative permeability of phase j, and s∗j is the normalized saturation of phase j, defined as

s∗j =

s j − sr j , n p 1 − j=1 sr j

(59)

where srj is the residual saturation of phase j. The parameters used to calculate the relative permeability curves are given in Table 1, and the curves are presented in Fig. 1. Water viscosity containing carbon dioxide was calculated using the Vogel-Fulcher-Tammann (VFT) correlation (Angell et al., 20 0 0). Friction theory (f-theory), which relates the viscosity of a fluid to an equation of state (Quiñones-Cisneros et al., 20 0 0; 20 01), was used to model the viscosity of the oil phase, which is composed by n-hexadecane and dissolved carbon dioxide. All other oil properties were calculated considering it as pure n-hexadecane.

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

187

Four different solution structures were found, which are defined by the ratio between the values of Henry’s constant of carbon dioxide in oil phase and Henry’s constant of carbon dioxide in water phase (KD ). If this ratio is smaller than 0.19, we call it Case I. When the value of KD is found in the range 0.19 < KD ≤ 0.457 it is defined as Case II, and if 0.457 < KD < 1 it is defined as Case III. If KD is greater than 1 (the last option), the solution is classified as Case IV. It is important to mention that Harveys´ correlation was used to determine Henrys´ constant only for Case I, for the other cases the relation KD was arbitrarily defined in order to find out different solution structures. 3.1. Case I: KD < 0.19 In addition to the data given in Table 2, for this case the ratio between the values of calculated Henry’s constants is KD = 0.02456. The structure of the solution (Figs. 2 and 3) is given by J − a → b − c → d − e → I, where (− ) denotes a rarefaction wave and ( → ) indicates a shock wave. The solution starts at point (J) with the first family rarefaction wave (concentration and temperature are constant), where oil sat-

Fig. 1. Relative permeability curves. Table 1 Relative permeability parameters. Property

Oil

Water

srj k0r j nj

0 1 2

0 1 2

From Eq. (9) and using Darcy’s law (Eq. (10)), we get the following expression to calculate the fractional flow of phase j:

kr j

fj =

μj n p k r j j=1

,

(60)

μj

where the relative permeabilities are calculated by Eq. (58) and the oil phase viscosity by f-theory model and water phase viscosity by the VFT correlation. Other parameters used to build the solution are given in Table 2. Specific heat at constant pressure data were obtained in Poling et al. (2001) and densities of pure components were calculated with Peng-Robinson equation of state (Robinson and Peng, 1978). Heat capacity and density of the rock are found in Schön (2011). For all cases cˆso = 0.0 0 02 volume fraction was set as initial solvent concentration. Fig. 2. Solution path: plane (so , fo ) for Case I. Table 2 Physical properties. Property

Symbol

Value

Unit

Reservoir pressure Injection pressure Water critical pressure Water critical temperature Oil saturation at injection condition Oil saturation at initial condition Solvent concentration at injection condition Solvent concentration at initial condition Temperature at injection condition Temperature at initial condition Oil heat capacity Water heat capacity Solvent heat capacity Rock heat capacity

P P (J) Pc Tc J so( ) I so( ) (J ) cˆso (I ) cˆso T (J) T (I) Cpo Cpw Cps Cpr

30.0 35.0 22.064 647.14 0.0 1.0 0.01 0.0 0 02 403.15 353.15 1459.37 4179.30 615.70 720.0

MPa MPa MPa K

Oil density Water density Solvent density Rock density Porosity

ρo ρw ρs ρr φ

628.22 803.27 657.15 2200.0 0.3

m3 m3 m3 m3 m3 m3 m3 m3

K K J KgK J KgK J KgK J kgK kg m3 kg m3 kg m3 kg m3 m3 m3

Fig. 3. Solution path: physical plane (xD , tD ) for Case I.

188

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

Fig. 4. Profiles: oil saturation, concentration and normalized temperature for Case I.

Fig. 5. Zoom in the shock trajectory (a → b) for Case I. Fig. 7. Phase space for so , cˆso and T for Case I.

() ditions fractional flow to T(J) and cˆso fractional flow at point (b). Next, there is another constant state zone (Fig. 3) followed by an oil saturation rarefaction (first family - concentration and temperature are constant) up to point (c). From (c) there is an oil saturation and temperature shock (third family - concentration remains constant) connecting intermediate fractional flow function to initial conditions fractional flow at point (d), the degenerate shock wave. This shock is followed by a constant state region (Fig. 3). After this constant state the last oil saturation rarefaction starts (first family concentration and temperature stay constant), and goes up to point (e), which is connected to initial conditions (I) through a Buckley– Leverett type shock (Buckley and Leverett, 1942). The solution profiles are presented in Fig. 4, where TD is the dimensionless temperature, TD = T /T (J ) . The shocks (a → b), the saturation-concentration shock, and (c → d), the degenerate shock, are zoomed in Figs. 5 and 6, respectively. In Fig. 7 we present the solution path in phase space for so , cˆso and T, where the continuous line represents the rarefaction waves and the dotted line denotes the shock waves. Segment J − a is the oil saturation rarefaction, a → b is the oil saturation and concentration shock, path b − c is another oil saturation rarefaction, c → d is the degenerate shock wave, segment d − e is the last oil saturation rarefaction, and from e → I is the Buckley– Leverett type shock (oil saturation shock). In Table 3 we can see I

Fig. 6. Zoom in the shock trajectory (c → d) for Case I.

uration changes from injection condition up to (a), where there is a constant state zone. At (a) there is an oil saturation and concentration shock (second family - temperature is constant while the oil saturation and concentration change) connecting injection con-

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

189

Table 3 Solution coordinates and dimensionless velocities for Case I.

J a b c d e I

Saturation so

Concentration cˆso

Temperature T

Dimensionless velocity

0.0 0 0 0 0.1938 0.2448 0.3509 0.4383 0.5927 1.0 0 0 0

0.0100 0.0100 0.0 0 02 0.0 0 02 0.0 0 02 0.0 0 02 0.0 0 02

403.15 403.15 403.15 403.15 353.15 353.15 353.15

0.0 0 0 0 0.1787 − 0.1796 0.1796 − 0.2485 0.5170 0.5170 − 0.7825 1.7276 1.7276

λ

Fig. 8. Solution path: plane (so , fo ) for Case II.

Fig. 11. Zoom in the shock trajectory (a → b) for Case II.

the coordinates of the solution and their respective dimensionless velocities. 3.2. Case II: 0.19 < KD ≤ 0.457

Fig. 9. Solution path: physical plane (xD , tD ) for Case II.

Besides the data given in Table 2, in this case the ratio between the values of Henry’s constants is set to KD = 0.2456. The structure of the solution (Figs. 8 and 9) is: J − a → b → c − d → I. The solution begins at point (J) with an oil saturation rarefaction (first family - concentration and temperature stay constant) up to point (a), followed by a third family shock wave

Fig. 10. Profiles: oil saturation, concentration and normalized temperature for Case II.

190

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195 Table 4 Solution coordinates and dimensionless velocities for Case II.

J a b c d I

Saturation so

Concentration cˆso

Temperature T

Dimensionless velocity

0.0 0 0 0 0.3380 0.4199 0.4653 0.5927 1.0 0 0 0

0.0100 0.0100 0.0100 0.0 0 02 0.0 0 02 0.0 0 02

403.15 403.15 353.15 353.15 353.15 353.15

0.0 0 0 0 0.5144 0.5144 − 0.5755 0.5755 − 0.9147 1.7274 1.7274

λ

Fig. 12. Zoom in the shock trajectory (b → c) for Case II.

Fig. 14. Solution path: plane (so , fo ) for Case III.

Fig. 13. Phase space for so , cˆso and T for Case II.

(degenerate shock - constant concentration), where oil saturation and temperature change. This shock connects injection conditions fractional flow function (point (a)) to the intermediate fractional (J ) flow function with parameters T(I) and cˆso at point (b). Next, there is a constant state zone (Fig. 9). Then, there is an oil saturation and concentration shock (second family - temperature remains constant) connecting point (b) to point (c) at initial conditions fractional flow curve, where a second constant state region begins (Fig. 9). At this point, there is another oil saturation rarefaction (concentration and temperature are constant) up to (d), which is connected to initial conditions (I) through a Buckley–Leverett type shock (Buckley and Leverett, 1942). The solution profiles are shown in Fig. 10. In Figs. 11 and 12 there is a zoom in shock regions (a → b), the degenerate shock, and (b → c), the saturationconcentration shock, respectively. We show the phase space solution path for this case in Fig. 13. From point J to point a there is an oil saturation rarefaction, segment a → b is the oil saturation and temperature shock (degenerate wave), b → c is the oil saturation and concentration shock, c − d is the second oil saturation rarefaction, and from d to I we have the Buckley–Leverett type shock. The coordinates of the solution and the respective dimensionless velocities can be found in Table 4.

Fig. 15. Solution path: physical plane (xD , tD ) for Case III.

3.3. Case III: 0.457 < KD < 1.0 In this case the ratio between the values of Henry’s constants is KD = 0.6. The structure of the solution path (Figs. 14 and 15) is: J − a → b − c → d − e → I. The solution starts at point (J) with the first family rarefaction (concentration and temperature are constant) wave, where oil saturation changes from injection condition up to point (a). From (a), there is an oil saturation and temperature

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

191

Fig. 16. Profiles: oil saturation, concentration and normalized temperature for Case III.

Fig. 18. Phase space for so , cˆso and T for Case III.

Fig. 17. Zoom in the fractional flow curves in regions (a → b − c ) for Case III.

shock (third family - concentration remains constant) connecting point (a) to the intermediate fractional flow curve with parame(J ) ters T(I) and cˆso at point (b). This shock is followed by a constant state zone (Fig. 15). An oil saturation rarefaction (first family - concentration and temperature are constant) wave begins at (b) up to point (c), where there is another constant state region (Fig. 15). Next, an oil saturation and concentration shock (second family temperature is constant) connects point (c) to point (d) at initial conditions fractional flow curve. This shock is followed by a constant state region (Fig. 15), and after this region there is another oil saturation rarefaction (concentration and temperature remain constant) up to point (e), and finally a Buckley–Leverett type shock connecting (e) to initial conditions point (I). The solution profiles are given in Fig. 16. Fig. 17 zooms region (a → b − c ), the degenerate shock and the following oil saturation rarefaction. In Fig. 18 the phase space solution path for this case is presented. Similarly to the previous cases, segment J − a is an oil saturation rarefaction, a → b is an oil saturation and temperature shock, b − c is the second oil saturation rarefaction, c → d is the oil saturation and concentration shock wave, d − e is the last oil saturation rarefaction, and from e → I we find the Buckley–Leverett type shock. In Table 5 we show the coordinates of the solution and the respective dimensionless velocities for this case.

Table 5 Solution coordinates and dimensionless velocities for Case III.

J a b c d e I

Saturation so

Concentration cˆso

Temperature T

Dimensionless velocity

0.0 0 0 0 0.3359 0.4110 0.4318 0.5264 0.5925 1.0 0 0 0

0.0100 0.0100 0.0100 0.0100 0.0 0 02 0.0 0 02 0.0 0 02

403.15 403.15 353.15 353.15 353.15 353.15 353.15

0.0 0 0 0 0.5131 0.5131 − 0.7403 0.8352 − 0.8364 0.8364 − 1.2708 1.7271 1.7271

λ

3.4. Case IV: KD > 1 The ratio between the values of Henry’s constants for this case is KD = 2.456. The structure of the solution path (Figs. 19 and 20) for this case is given by: J − a → b − c − d − e → I. Therefore, the solution starts at injection conditions point (J) with an oil saturation rarefaction of the first family (concentration and temperature remain constant) up to point (a). Next, an oil saturation and temperature shock wave, the degenerate shock (third family - concentration is constant), connects the injection to intermediate fractional flow (J ) functions with parameters T(I) and cˆso at point (b), followed by a

192

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

Fig. 22. Zoom in the shock trajectory (a → b) for Case IV.

Fig. 19. Solution path: plane (so , fo ) for Case IV.

constant state region (Fig. 20). From (b), there is another oil saturation rarefaction (first family - concentration and temperature are constant) up to point (c). Then, an oil saturation and concentration rarefaction wave (second family - temperature remains constant) connects point (c) to point (d) (initial conditions fractional flow function). This rarefaction is followed by a constant state region (Fig. 20). From (d) there is the last oil saturation rarefaction (concentration and temperature are constant) up to (e), and a Buckley– Leverett type shock linking this point to initial conditions (I). In Fig. 21 the solution profiles are presented. In Fig. 22 the degenerate shock (a → b) trajectory is presented in detail. The phase space solution path for so , cˆso and T is shown in Fig. 23. Again, J − a is an oil saturation rarefaction, a → b is the degenerate shock wave, b − c is the second oil saturation rarefaction, c − d is the oil saturation and concentration rarefaction wave, d − e is the last oil saturation rarefaction, and e → I is the oil saturation shock (Buckley–Leverett). The coordinates of the solution and the respective dimensionless velocities for this case are presented in Table 6. Discussion Fig. 20. Solution path: physical plane (xD , tD ) for Case IV.

The main difference between the solutions of the four cases is the relation between the temperature and the concentration wave speeds. For cases II, III and IV, the temperature shock speed

Fig. 21. Profiles: oil saturation, concentration and normalized temperature for Case IV.

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

193

Fig. 25. Recovery factor for carbonated water injection at cˆsw = 0.004 for different KD compared to hot waterflooding (HW). Table 7 Physical properties for sensitivity analysis.

Fig. 23. Phase space for so , cˆso and T for Case IV. Table 6 Solution coordinates and dimensionless velocities for Case IV.

J a b c d e I

Saturation so

Concentration cˆso

Temperature T

Dimensionless velocity

0.0 0 0 0 0.3249 0.3538 0.4873 0.5553 0.5920 1.0 0 0 0

0.0100 0.0100 0.0100 0.0100 0.0 0 02 0.0 0 02 0.0 0 02

403.15 403.15 353.15 353.15 353.15 353.15 353.15

0.0 0 0 0 0.5064 0.5064 − 0.5971 1.2486 1.2810 − 1.4668 1.7255 1.7255

λ

(degenerate shock) is slower than the concentration wave speed, while for Case I the concentration wave (shock) is slower. For this reason, most part of the solution for Case I takes place at injection temperature (higher temperature), which leads to lower carbon dioxide solubility in the oil phase (Nasehi and Asghari, 2010). In Case II, there is a constant state before the concentration shock, and in Case III there is a saturation rarefaction before the concentration shock. Case IV presents a saturation rarefaction followed by a concentration rarefaction. The remaining part of the solution is the same for these cases (II, III and IV). The oil production rate for all the discussed cases is presented in Fig. 24 and the recovery factor (RF) of the solutions is compared to hot waterflooding in Fig. 25. According to Fig. 25, intermediate values of KD present the best results. High Henry’s constants ratios indicate small CO2 transfer from carbonated water to oil, leading to a similar (or even slightly worse) recovery factor compared to hot waterflooding. On the other hand, very small ratios turn the injec-

Fig. 24. Oil production rate for carbonated water injection for different KD .

Property

Symbol

Value

Unit

Reservoir pressure Injection pressure Mean pressure Water critical pressure Water critical temperature Mean temperature Oil saturation at injection condition Oil saturation at initial condition Solvent concentration at initial condition Temperature at injection condition Temperature at initial condition Oil heat capacity Water heat capacity Solvent heat capacity Rock heat capacity

P P (J) P¯ Pc Tc T¯

30.0 35.0 32.5 22.064 647.14 378.15 0.0 1.0 0.0 0 02 403.15 353.15 1459.37 4179.30 615.70 720.0

MPa MPa MPa MPa K K

Oil density Water density Solvent density Rock density Porosity Injection rate Reservoir cross sectional area Reservoir length Injection time

so( ) I so( ) (I ) cˆso T (J) T (I) Cpo Cpw Cps Cpr J

ρo ρw ρs ρr φ q A L tf

628.22 803.27 657.15 2200.0 0.3 50.0 10 0.0 0 10 0 0.0 0 3

m3 m3 m3 m3 m3 m3

K K J KgK J KgK J KgK J kgK kg m3 kg m3 kg m3 kg m3 m3 m3 m3 d 2

m m years

tion carbonated water into hot water, because the carbon dioxide is mostly dissolved in oil. When most of the CO2 is transferred to the oil, it causes reduction in oil viscosity as well as a larger oil saturation due to oil swelling. So, we can conclude that the mass transfer of CO2 plays an important role in this EOR process, due to oil swelling and viscosity reduction. 4. Sensitivity analysis In this section, we evaluate the effect of different carbon dioxide concentrations in injected water. All phases properties and rock-fluid interactions were calculated using the same procedure as described in Section 3. Reservoir and operational characteristics are found in Table 7. The initial solvent concentration was set to cˆso = 0.0 0 02 volume, and the calculated ratio between Henry’ s constants was KD = 0.02456 (Case I). Three different carbon dioxide concentrations at injection conditions were analyzed: cˆsw = 0.001, cˆsw = 0.004, and cˆsw = 0.007. The structure of the solution path for these injection conditions is the same as the Case I presented previously (Figs. 2 and 3). For the smallest and greatest concentrations, the coordinates of the solution and the respective dimensionless velocities can be found in Tables 8 and 9. The oil production rate for the three concentrations

194

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

Table 8 Solution coordinates and dimensionless velocities for carbonated water at 0.1% vol.

is shown in Fig. 26 and the RF can be found in Fig. 27. These results are also compared to hot waterflooding. According to Fig. 27, as the solvent concentration in the water phase increases, the production rate and the recovery factor also rise.

combined thermal-solvent EOR method. The governing system of equations (hyperbolic system) is composed by three equations and it is solved using the method of characteristics. The solution path is composed by rarefaction and shock waves and constant states, and this problem presents a degenerate wave. Along one of the waves oil saturation changes while the solvent concentration in the oil phase and the temperature are constant. Oil saturation and concentration changes along another wave, while the temperature is constant. For the last wave, the oil saturation and temperature change while the concentration remains constant. The Rankine–Hugoniot conditions were used to build up the shock conditions for the three different shock types. Four different solutions were presented for a particular set of reservoirfluid properties, the difference among them depend on the ratio between the value of Henry’s constant of the solvent in oil phase and the value of Henry’s constant of the solvent in water phase (KD ). In one of the cases (KD < 0.19), the solution is divided into seven regions, consisting of an oil saturation rarefaction (first family), followed by two constant state zones separated by an oil saturation and concentration shock (second family). Next, there is an oil saturation rarefaction and then a temperature shock appears (degenerate wave). This shock is followed by a third constant state region, an oil saturation rarefaction wave, and the first type shock (Buckley–Leverett type). When 0.19 < KD ≤ 0.457, the solution is divided into six regions, composed by an oil saturation rarefaction (first family), followed by an oil saturation and temperature shock, the degenerate wave (third family). Next, there are two constant state regions separated by an oil saturation and concentration shock (second family). Then, a second constant region is followed by another oil saturation rarefaction, and finally the Buckley–Leverett type shock. For another case (0.457 < KD < 1), the solution is divided into seven regions. It is composed by an oil saturation rarefaction (first family) followed by the degenerate wave (third family). Then, there are two constant state regions separated by an oil saturation rarefaction (first family). The second constant region is followed by an oil saturation and concentration shock (second family). Following this wave there is the last constant state zone, another oil saturation rarefaction, and finally the Buckley–Leverett type shock. If KD > 1, the solution is composed by seven regions. It consists of an oil saturation rarefaction (first family) followed by a temperature shock (degenerate wave) and a constant state region. Next, an oil saturation rarefaction wave appears before an oil saturation and concentration rarefaction wave (second family). This wave is followed by another constant state zone, then an oil saturation rarefaction, and the first type shock (Buckley–Leverett type). It is shown that the mass transfer of CO2 plays an important role in this EOR process, due to oil swelling and viscosity reduction. As the ratio KD increases the recovery factor (RF) also increases, because more solvent is transferred to oil phase. However, if KD exceeds 1, there is so much carbon dioxide transfer to the oil phase that this process becomes similar to hot waterflooding, leading to the smallest recovery factor. A sensitivity analysis regarding the carbon dioxide concentration showed that for these specific conditions the higher the concentration the higher the recovery factor. This solution may be used for fast and efficient evaluation of technical aspects of the best Enhanced Oil Recovery method for a particular field and for the development of streamline simulators for systems where diffusion effects are negligible.

5. Conclusions

Acknowledgements

In this work, we develop an analytical solution for the problem of one-dimensional oil displacement by hot carbonated water as a

Tamires Marotto would like to thank the Brazilian Government Agency CAPES for the support of her master’s degree. The authors

J a b c d e I

Saturation so

Concentration cˆso

Temperature T

Dimensionless velocity

0.0 0 0 0 0.1987 0.2219 0.3509 0.4383 0.5927 1.0 0 0 0

0.0026 0.0026 0.0 0 02 0.0 0 02 0.0 0 02 0.0 0 02 0.0 0 02

403.15 403.15 403.15 403.15 353.15 353.15 353.15

0.0 0 0 0 0.1756 − 0.1758 0.1758 − 0.2082 0.5170 0.5170 − 0.7825 1.7276 1.7276

λ

Table 9 Solution coordinates and dimensionless velocities for carbonated water at 0.7% vol.

J a b c d e I

Saturation so

Concentration cˆso

Temperature T

Dimensionless velocity

0.0 0 0 0 0.1880 0.2626 0.3509 0.4383 0.5927 1.0 0 0 0

0.0193 0.0193 0.0 0 02 0.0 0 02 0.0 0 02 0.0 0 02 0.0 0 02

403.15 403.15 403.15 403.15 353.15 353.15 353.15

0.0 0 0 0 0.1825 − 0.1843 0.1843 − 0.2838 0.5170 0.5170 − 0.7825 1.7276 1.7276

λ

Fig. 26. Oil production rate for carbonated water injection at different concentrations.

Fig. 27. Recovery factor for carbonated water injection at different concentrations compared to hot waterflooding (HW).

T.A. Marotto and A.P. Pires / International Journal of Multiphase Flow 115 (2019) 181–195

also thank Petrobras for funding this research through the SIGER Research Network under contract number 0 050 0 082421139. References Angell, C.A., Ngai, K.L., McKenna, G.B., McMillan, P.F., Martin, S.W., 20 0 0. Relaxation in glassforming liquids and amorphous solids. J. Appl. Phys. 88 (6), 3113– 3157. Bedrikovetsky, P., 1993. Mathematical Theory of Oil and Gas Recovery. Kluwer Academic Publishers, Londres. Bruining, J., Marchesin, D., 2006. Analysis of nitrogen and steam injection in a porous medium with water. Transp. Porous Media 62 (3), 251–281. doi:10.1007/ s11242- 005- 1733- 2. Buckley, S.E., Leverett, M.C., 1942. Mechanisms of fluid displacement in sands. SPE AIME 146 (1), 107–116. SPE-942107. Christensen, R.J., 1961. Carbonated waterflooding results - Texas and Oklahoma. Annual Meeting of Rocky Mountain Petroleum Engineers of AIME. Society of Petroleum Engineers, Farmington, New Mexico SPE-66. 10.2118/66-MS. Claridge, E.L., Bondor, P.L., 1974. A graphical method for calculating linear displacement with mass transfer and continuously changing mobilities. Soc. Petrol. Eng. J. 14 (06), 609–618 SPE-4673. 10.2118/4673-PA. Corey, A.T., Rathjens, C.H., Henderson, J.H., Wyllie, M.R.J., 1956. Three-phase relative permeability. J. Petrol. Technol. 8 (11), 63–65. SPE-737. Dindoruk, B., 1992. Analytical Theory of Multiphase Multicomponent Displacement in Porous Media. Stanford University Doctorate. Dindoruk, D.S., Dindoruk, B., 2008. Analytical and numerical solution of nonisothermal Buckley-Leverett flow including tracers. SPE Reserv. Eval. Eng. 11 (3), 555–564. SPE-102266. Dong, Y., Dindoruk, B., Ishizawa, C., Lewis, E., Kubicek, T., 2011. An experimental investigation of carbonated water flooding. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers, Denver, Colorado. SPE-145380. Ezekwe, N., 2011. Petroleum Reservoir Engineering Practice. Pearson Education, Boston. Fayers, F.J., Perrine, R.L., 1958. Mathematical description of detergent flooding in oil reservoirs. Fall Meeting of the Society of Petroleum Engineers of AIME. Society of Petroleum Engineers, Houston, Texas SPE-1132G. 10.2118/1132-G. Fayers, F.J., 1962. Some theoretical results concerning the displacement of a viscous oil by a hot fluid in a porous medium. J. Fluid Mech. 13 (5), 65–76. Green, D.W., Willhite, G.P., 1998. Enhanced Oil Recovery. Society of Petroleum Engineers, Richardson. Harvey, A.H., 1996. Semiempirical correlation for Henry‘s constants over large temperature ranges. AlChE J. 42 (5), 1491–1494. Hickok, C.W., Ramsey Jr., H.J., 1962. Case histories of carbonated waterflood in Dewey-Bartlesville fields. SPE Secondary Recovery Symposium. Society of Petroleum Engineers, Wichita Falls, Texas SPE-333. 10.2118/333-MS. Holke, D.C., Huebner, W.B., 1971. Thermal stimulation and mechanical techniques permit increased recovery from unconsolidated viscous oil reservoir. SPE Anadarko Basin/Oklahoma City Joint Meeting. Society of Petroleum Engineers, Oklahoma City, Oklahoma. SPE-3671. Kechut, N.I., Jamiolahmady, M., Sohrabi, M., 2011. Numerical simulation of experimental carbonated water injection (CWI) for improved oil recovery and CO2 storage. J. Petrol. Sci. Eng. 77 (1), 111–120. doi:10.1016/j.petrol.2011.02.012. Kokal, S., Al-Kaabi, A., 2010. Enhanced oil recovery: challenges & opportunities. World Petroleum Coun. 64 (1), 64–69. Lake, L.W., 1989. Enhanced Oil Recovery. Prentice-Hall, Englewood Cliffs, NJ. Lambert, W., Marchesin, D., Bruining, J., 2010. The Riemann solution for the injection of steam and nitrogen in a porous medium. Transp. Porous Media 81 (3), 505–526. doi:10.10 07/s11242-0 09-9419-9. LeVeque, R.J., 2002. Finite Volume Method for Hyperbolic Problems. Cambridge University Press, Cambridge. Martin, W.L., Dew, J.N., Powers, M.L., Steves, H.B., 1968. Results of a tertiary hot waterflood in a thin sand reservoir. J. Petrol. Technol. 20 (7), 739–750. SPE-1918. Mosavat, N., Torabi, F., 2014. Performance of secondary carbonated water injection in light oil systems. Ind. Eng. Chem. Res. 53 (3), 1262–1273.

195

Nasehi, M., Asghari, K., 2010. Use of CO2 in heavy-oil waterflooding. In: SPE International Conference on CO2 Capture, Storage, and Utilization. Society of Petroleum Engineers, New Orleans, Louisiana SPE-139672. 10.2118/139672-ms. Perez, J.M., Poston, S.W., Sharif, Q.J., 1992. Carbonated water imbibition flooding: an enhanced oil recovery process for fractured reservoirs. SPE/DOE Enhanced Oil Recovery Symposium. Society of Petroleum Engineers, Tulsa, Oklahoma SPE24164-MS. 10.2118/24164-MS. Picha, M.S., 2007. Enhanced oil recovery by hot CO2 flooding. In: SPE Middle East Oil & Gas Show and Conference. Society of Petroleum Engineers, Manama, Bahrain. SPE-105425 Poling, B.E., Prausnitz, J.M., O’Connell, J.P., 2001. The Properties of Gases and Liquids, 5th ed. McGraw-Hill. Polyanin, A.D., Zaitsev, V.F., 2003. Handbook of Nonlinear Partial Differential Equation. Chapman & Hall/CRC, Boca Raton. Pope, G.A., 1980. The application of fractional flow theory to enhanced oil recovery. Soc. Petrol. Eng. J. 20 (3), 191–205. SPE-7660. Prausnitz, J.M., Lichtenthaler, R.N., Azevedo, E.G., 1986. Molecular Thermodynamics of Fluid-Phase Equilibria. Prentice-Hall, Englewood Cliffs, NJ. Quiñones-Cisneros, S.E., Zéberg-Mikkelsen, C.K., Stenby, E.H., 20 0 0. The friction theory (f-theory) for viscosity modeling. Fluid Phase Equilib. 169 (2), 149–276. Quiñones-Cisneros, S.E., Zéberg-Mikkelsen, C.K., Stenby, E.H., 2001. The friction theory for viscosity modeling: extension to crude oil systems. Chem. Eng. Sci. 56 (24), 7007–7015. Riazi, M., Sohrabi, M., Jamiolahmady, M., Ireland, S., Brown, C., 2009. Oil recovery improvement using CO2 -enriched water injection. In: EUROPEC/EAGE Conference and Exhibition. Society of Petroleum Engineers, Amsterdam, The Netherlands. SPE-121170-MS. Riazi, M., Jamiolahmady, M., Sohrabi, M., 2011a. Theoretical investigation of pore-scale mechanisms of carbonated water injection. J. Petrol. Sci. Eng. 75 (3–4), 312–326. Riazi, M., Sohrabi, M., Jamiolahmady, M., 2011b. Experimental study of pore-scale mechanisms of carbonated water injection. Transp. Porous Media 86 (1), 73–86. Robinson, D.B., Peng, D.Y., 1978. The Characterization of the Heptanes and Heavier Fractions for the GPA Peng-Robinson Programs, Research Report RR-28 ed. Gas Processors Association, Tulsa, Okla. Saul, A., Wagner, W., 1987. International equations for the saturation properties of the ordinary water substance. J. Phys. Chem. 16 (4), 893–901. Schön, J., 2011. Physical Properties of Rocks, vol. 8, first ed. Elsevier, Amsterdam, The Netherlands. Seto, C., 2007. Analytical Theory for Two-Phase, Multicomponent Flow in Porous Media with Adsorption. Stanford University Doctorate. Shakiba, M., Ayatollahi, S., Riazi, M., 2016. Investigation of oil recovery and CO2 storage during secondary and tertiary injection of carbonated water in an iranian carbonate oil reservoir. J. Petrol. Sci. Eng. 137 (1), 134–143. doi:10.1016/j.petrol. 2015.11.020. Sohrabi, M., Riazi, M., Jamiolahmady, M., Kechut, N.I., Ireland, S., Robertson, G., 2011. Carbonated water injection (CWI) - a productive way of using CO2 for oil recovery and CO2 storage. Energy Procedia 4 (1), 2192–2199. Sohrabi, M., Kechut, N.I., Riazi, M., Jamiolahmady, M., Ireland, S., Robertson, G., 2012. Coreflooding studies to investigate the potential of carbonated water injection as an injection strategy for improved oil recovery and CO2 storage. Transp. Porous Media 91 (1), 101–121. Spillette, A.G., Nielsen, R.L., 1968. Two-dimensional method for predicting hot waterflood recovery behavior. J. Petrol. Technol. 20 (6), 627–638. SPE-1895. Spillette, A.G., 1965. Heat transfer during hot fluid injection into an oil reservoir. J. Can. Pet. Technol. 4 (4), 213–218. PETSOC-65-04-06. Welge, H.J., 1952. A simplified methods for computing oil recovery by gas or water drive. J. Petrol. Technol. 4 (04), 91–98 SPE-124. 10.2118/124-G. Willman, B.T., Valleroy, V.V., Runberg, G.W., Cornelius, A.J., Powers, L.W., 1961. Laboratory studies of oil recovery by steam injection. J. Petrol. Technol. 13 (7), 681–690. SPE-1537. Zhao, L., 2004. Steam alternating solvent process. SPE International Thermal Operations and Heavy Oil Symposium and Western Meeting. Society of Petroleum Engineers, Bakersfield, California. SPE-86957. Zhu, J., 2003. Multicomponent Multiphase Flow in Porous Media with Temperature Variation or Adsorption. Stanford University Doctorate.