Modeling Transient Flow in CO₂ Injection Wells by Considering the Phase Change

Abstract

A transient model to simulate the temperature and pressure in CO₂ injection wells is proposed and solved using the finite difference method. The model couples the variability of CO₂ properties and conservation laws. The maximum error between the simulated and measured results is 5.04%. The case study shows that the phase state is primarily controlled by the wellbore temperature. Increasing the injection temperature or decreasing the injection rate contributes to obtaining the supercritical state. The variability of density can be ignored when the injection rate is low, but for a high injection rate, ignoring this may cause considerable errors in pressure profiles.

1. Introduction

CO₂ has been widely utilized for oil and gas exploitation. It can be injected into the reservoir to enhance oil recovery (EOR), as it can sweep the mobilized oil to the surface [1,2]. Additionally, in some unconventional reservoirs that require fracturing work, CO₂ is adopted as a fracturing fluid since it causes little damage to the reservoir [3,4]. In most cases, the injected CO₂ is in a supercritical state. When temperature and pressure exceed the critical value (31.16 °C; 7.382 MPa), the CO₂ viscosity and interface tension are similar to those of gas, but the density is similar to that of liquid, which is helpful for effective oil flooding or reservoir fracturing. Thus, the phase state of the CO₂ in injection wells has attracted many researchers’ attention in recent years [5,6]. The key to controlling the CO₂ phase is to simulate the temperature and pressure of flowing CO₂ in the wellbore.

On account of the complex well-bottom conditions, the prediction of CO₂ temperature and pressure is a challenge. Most pressure simulation works have focused on the friction coefficient of CO₂ and the friction pressure drop. Li et al. [7] and Wang et al. [8] studied the pressure drop in CO₂ fracturing wells using numerical methods and experiments. They proposed empirical relations to calculate the CO₂ friction pressure drop in CO₂ wells. For temperature simulation, most temperature simulation models were based on the following approaches: the semi-steady state method and transient model. The semi-steady state model assumes a steady flow in the wellbore but transient heat conduction in the formation [9,10], while the transient model regards the fluid flow in the wellbore as a transient process [11,12,13]. Compared with the semi-steady model, the transient model is more complicated to use but can accurately simulate the unsteady flow. In 1982, Cranshaw and Bolling [14] performed a numerical simulation of non-isothermal flow in CO₂ wells by considering two-phase flow. Zhang and Tang [15] adopted the semi-steady model to predict temperature and pressure in CO₂ injection wells. Field data with a low injection rate were used to verify the model. Yasunami [16] proposed a numerical method to simulate the CO₂ flow and suggested controlling the CO₂ phase by adjusting the injection rate. Paterson et al. [17] adopted Hasan’s model to estimate the wellbore temperature. Lu and Connell [18] proposed a transient model to study CO₂ flow during geological storage. In the last 5 years, many scholars have advanced the prediction model from different aspects. Sun et al. [19] and Song et al. [20] discussed the CO₂ temperature in a wellbore by using the semi-steady state model. The Span–Wagner state equation was used in their models to calculate CO₂ physical properties. Yi et al. [21], Guo and Zeng [22], Gong et al. [23], Wang et al. [24], Lyu et al. [25,26], Li et al. [27], and Yang et al. [28] proposed transient models to simulate CO₂ temperature and discussed the sensitivity of some key factors.

To improve the understanding of CO₂ flow in injection wells, the primary goal of this work was to study the phase state in CO₂ injection wells. A transient model considering the CO₂ phase transition was established by the conservation law. The variability of physical properties caused by the phase change, mass conservation, momentum conservation, and energy conservation was coupled to obtain the wellbore temperature and pressure. The control methods of phase state were investigated using case studies. Finally, the influences of density variation caused by the phase change on the flow were analyzed.

2. Transient Method

2.1. Establishment of the Model

Figure 1 is a schematic of a CO₂ injection well. In this well, CO₂ is injected by tubing. In some cases, CO₂ can also be injected by the annulus between tubing and production casing [23]. Overall, the following work is based on the tubing injection. The simulation model of annular injection can be developed by the same method presented below. To simplify the model, the following assumptions were made: (a) the axial heat conduction is ignored [23]; (b) the wellbore is full of flowing CO₂ before injection, and the initial CO₂ has attained heat balance with surrounding formation before CO₂ injection; (c) the fluid velocity, temperature, and pressure are constant at the same depth; (d) the flow and the well depth are one-dimensional. Basic equations were established by the conservation laws, as presented below. The derivation of the following relations is presented in Appendix A.

• Mass conservation

The continuity equation derived by mass conservation is shown in Equation (1).

$$\frac{\partial \rho }{\partial t}+V\frac{\partial \rho }{\partial z}+\rho \frac{\partial V}{\partial z}=0$$

(1)

• Momentum conservation

As the model is assumed to be one-dimensional, forces on the liquid consist of the gravity, the pressure, and the friction caused by the viscosity. By the momentum conservation, the pressure gradient can be expressed as Equation (2). The first term on the right side indicates the gravity pressure drop, the second term accounts for the acceleration pressure drop, and the third one is the friction pressure drop. Considering that the well trajectory may be curving or horizontal, the gravity acceleration ($g_z$) should be calculated by $g_z=g\cos \theta$.

$$\frac{\partial P}{\partial z}=\rho g_z-\rho \frac{\partial V}{\partial z}-\rho F_f$$

(2)

• Energy conservation

According to energy conservation, the fluid temperature can be described by Equation (3). The first term on the right side is the contribution of heat conduction. The second one is the heat caused by friction, and the third term denotes the influence of pressure on temperature.

$$\frac{c_{P}dT}{dt}=Q_c+Q_f+\left(\frac{1}{\rho }+c_j{c}_P\right)\frac{dP}{dt}$$

(3)

In Equations (1)–(3), V denotes the fluid velocity (m/s); z is the measured depth of the wellbore (m); $g_z$ is the gravity acceleration in the well trajectory (m/s²); $g$ represents the gravity acceleration (9.81 m/s²); $\theta$ is the deviation angle (rad); P is the fluid pressure (MPa); $\rho$ is the fluid density (Kg/m³); t is the injection time (s); $F_f$ is the friction per unit mass (N/Kg); $c_P$ is the specific heat (J/(K × Kg)); T is the fluid temperature (K); $Q_c$ represents the heat transfer by conduction (J/(s × Kg)); $Q_f$ is the heat generation rate due to friction per unit mass (J/(s × Kg)); $c_j$ is the Joule–Thomson coefficient (K/MPa).

2.2. Calculation of Key Parameters

2.2.1. Properties of CO₂

By considering the phase change, the physical properties of CO₂ are variable with temperature and pressure, including the density, viscosity, and heat capacity. REFPROP, a commercial piece of software developed by the National Institute of Standards and Technology (NIST), provides a convenient way to output accurate CO₂ properties. Therefore, this work adopted REFPROP to obtain CO₂ properties at different temperatures and pressures.

2.2.2. Heat Transfer

The wellbore can be regarded as a multiple-layer cylinder with no internal heat source. Inside the cylinder is flowing fluid, and outside the cylinder is formation. The convection heat transfer in Equation (3) can be calculated by Equation (4).

$$Q_c=\frac{(T_g-T)}{m{\displaystyle \sum R_t}}$$

(4)

$R_t$ is the thermal resistance of one layer, as shown in Equation (5).

$$R_t=\frac{1}{2\pi \lambda H}\ln \frac{r_0}{r_i}$$

(5)

In Equations (4) and (5), $T_g$ is the formation temperature (K); m denotes the fluid mass (Kg); $R_t$ is the thermal resistance of one layer (K/W); λ is the thermal conductivity (J/(m∙s∙K)); H is the height of fluid unit (m); $r_0$ is the outer radius of one layer (m); $r_i$ is the inner radius of one layer (m).

2.2.3. Heat Transfer

The friction of unit mass can be obtained by Equation (6).

$$F_f=\frac{f{V}^2}{4r_t}$$

(6)

$r_t$ is the inner radius of tubing (m); $f$ is the friction factor, which is related to the Reynold number. Chen’s relation was used to calculate the friction factor [28], as presented in Equation (7). This equation covers most Reynold numbers and has been widely used in CO₂ wells. In Equation (7), ∆ denotes the roughness of tubing (m); Re represents the Reynold number.

$$f={\left(-2\lg \left(\frac{\Delta /(2r_t)}{3.7065}-\frac{5.0452}{\mathrm{Re}}\lg \left(\frac{{(\Delta /(2r_t))}^{1.1098}}{2.8257}+\frac{5.8506}{{\mathrm{Re}}^{0.8981}}\right)\right)\right)}^{-2}$$

(7)

The heat caused by friction is the work of friction, as shown in Equation (8).

$$Q_f=F_{f}V$$

(8)

2.3. Initial and Boundary Conditions

2.3.1. Initial Conditions

The main initial parameters include the wellbore velocity profile, pressure profile, and temperature profile. The velocity, pressure, and temperature profiles were determined by the initial motion state of wellbore fluid. Assuming the initial CO₂ is flowing steadily, the initial conditions should be Equation (9).

$$\{\begin{array}{l}{V|}_{t=0}=V_S\\ {T|}_{t=0}=T_S\\ {P|}_{t=0}=P_S\end{array}$$

(9)

$V_S$ is the steady fluid velocity profile (m/s); $T_S$ is the steady fluid temperature profile (K); $P_S$ is the steady fluid pressure profile (MPa).

The steady flow profiles in Equation (9) can be calculated by neglecting the partial differential of time, as presented in Equations (10)–(12).

$$V\frac{\partial \rho }{\partial z}+\rho \frac{\partial V}{\partial z}=0$$

(10)

$$\frac{\partial P}{\partial z}=\rho g-\rho F_f$$

(11)

$$V\frac{c_{P}dT}{dz}=Q_c+Q_f+\left(\frac{1}{\rho }+c_j{c}_P\right)V\frac{dP}{dz}$$

(12)

2.3.2. Boundary Conditions

The boundary conditions are presented in Equations (13) and (14). The formation temperature was assumed constant outside the wellbore. Additionally, the injection temperature, pressure, and rate were constant at the wellhead. $V_{in}$ is the injection velocity (m/s); $P_{in}$ is the injection pressure (MPa); $T_{in}$ is the injection temperature (K).

$${T|}_{r=\infty }=T_g$$

(13)

$$\{\begin{array}{l}{V|}_{z=0}=V_{in}\\ {T|}_{z=0}=T_{in}\\ {P|}_{z=0}=P_{in}\end{array}$$

(14)

2.4. Solution Method

This work used the explicit method to solve the formulations considering its simplicity and short calculation time. Equations (15)–(17) are the explicit difference equations of the governing equations.

$$V_{i+1}^n=V_i^n-\frac{1}{{\rho }_i^n}\left(\frac{{\rho }_i^n-{\rho }_i^{n-1}}{\Delta t}+V_i^n\frac{{\rho }_{i+1}^n-{\rho }_i^n}{\Delta z}\right)\Delta z$$

(15)

$$P_{i+1}^n=P_i^n+P_i^n\Delta z\left(g{z}_i^n-\frac{V_i^n-V_i^{n-1}}{\Delta t}-V_i^n\frac{V_{i+1}^n-V_i^n}{\Delta z}-F_f{}_i^n\right)$$

(16)

$$\begin{array}{l}{T}_i^{n+1}=T_i^n+\Delta t\{\frac{Q_c{}_i^n+Q_f{}_i^n}{c_P{}_i^n}+\left(\frac{1}{{\rho }_i^n{c}_P{}_i^n}+c_j{}_i^n\right)\left(\frac{P_i^{n+1}-P_i^n}{\Delta t}+V_i^n\left(\frac{P_i^n-P_{i-1}^n}{\Delta z}\right)\right)\\ -\frac{V_i^n}{\Delta z}(T_i^n-T_{i-1}^n)\}\end{array}$$

(17)

The mesh was divided axially along with the well depth. As the difference equations are explicit, the time step should be controlled strictly to avoid divergence. The time step in this work was set to the fluid flow limit in one grid at each step. The calculation procedure is presented in Figure 2, where $T_a$ is the assumed temperature (K); $P_a$ is the assumed pressure (MPa); Vt is the trial velocity (m/s); $T_t$ is the trial temperature (K); $P_t$ is the trial pressure (MPa). In the following calculation, the space step was set to 50 m, and the time step was set to approximately 5 min. Using the CPU i7-8700, the following cases presented in this work can be performed with a calculation time of shorter than 5 s.

3. Validation of the Model

The field data in Cao-8 well [15] were used to validate the proposed model. In this well, the injection rate was relatively small as CO₂ was injected to displace oil in the formation. Two CO₂ injection tests were conducted in this well. The first test took 13 h, and the injection rate was 55.4 t/d. The injection pressure and temperature were 24.5 MPa and −20 °C, respectively. The measured field data and simulated data at the well bottom are listed in

Table 1. The measured and simulated data of Cao-8 in the first test.

Temperature			Pressure
Measured	Simulated	Relative	Measured	Simulated	Relative
data (°C)	data (°C)	errors (%)	data (MPa)	data (MPa)	errors (%)
101.00	101.90	0.89	52.02	52.51	0.94

. This shows that the maximum relative error of temperature is 0.94%, which is acceptable in engineering.

In the second injection test, the injection temperature and pressure were 20 °C and 30 MPa, respectively. The measured and simulated temperature and pressure in the tubing after 26.52 days of injection at the rate of 21.17 t/d are presented in Figure 3. Figure 3 shows that the discrepancies between simulated data and measured data are small. The error analysis shows that the maximum relative error is 5.04%. Thus, by comparing the field data from the Cao-8 well, the proposed model is verified.

4. Results and Discussion

4.1. Case Study—A Fracturing Well

The W-16 well is a typical fracturing well in the Jianghan Basin, where CO₂ is injected into the wellbore at the rate of 4 m³/s. The wellbore architecture is presented in Figure 1, and the specific parameters are listed in

Table 2. Parameters of W-16.

Parameters	Value	Parameters	Value
Measured depth (m)	2800	Wellbore diameter (mm)	444
Tubing external diameter (mm)	89	Tubing wall (mm)	6.5
Surface casing external diameter (mm)	339.7	Surface casing wall (mm)	9.65
Surface casing depth (m)	63	Cement outside surface casing	Surface
Technical casing external diameter (mm)	244.47	Technical casing wall (mm)	10.03
Technical casing depth (m)	969	Cement outside technical casing	Surface
Production casing external diameter (mm)	139.7	Production casing wall (mm)	9.17
Production casing depth (m)	2775	Cement outside production casing (m)	1424
Geothermal gradient (K/100m)	2.95	Surface temperature (K)	24
Injection rate (m3/min)	4	Injection pressure (MPa)	80
Injection temperature (K)	−20	Injection time (hour)	1
String thermal conductivity (J/(m∙s∙K))	53	Annular fluid thermal conductivity (J/(m∙s∙K))	0.557
Cement thermal conductivity (J/(m∙s∙K))	0.627	Formation thermal conductivity (J/(m∙s∙K))	1.6

. According to the transient model developed in Section 2, the pressure and temperature distributions are presented in Figure 4 and Figure 5.

Figure 4 shows that the pressure decreased significantly with the increase in well depth, which was caused by a considerable drop in the fraction pressure. However, with an increase in injection time, the pressure changed slightly. During the injection period, the wellbore pressure profiles at any time were greater than the supercritical pressure. These high-pressure profiles mean that the phase state of CO₂ in the tubing is either a liquid state or supercritical state, according to the phase state map of CO₂ [25]. Then, the phase state of CO₂ must be determined by the temperature distributions.

Figure 5 illuminates that the CO₂ temperature increased with the increase in well depth. If the temperature was more than the critical value, the CO₂ was the supercritical state. On the contrary, the temperature lower than the critical value denoted the liquid state, as presented in Figure 5. When the injection time was zero, the temperature profile was calculated by the initial conditions. This showed that CO₂ below 900 m was the supercritical state at the beginning. However, with the injection of low-temperature CO₂, the temperature profiles decreased with the injection time. This was caused by the heat exchange with the injected low-temperature CO₂. When CO₂ flowed steadily under the boundary conditions, the bottom temperature was below the critical temperature, so the CO₂ in the tubing was the liquid state.

In sum, in the fracturing well, the phase state of CO₂ was primarily attributed to the temperature distribution. This case study showed that the CO₂ in W-16 could not attain the supercritical state when the flow was steady. Thus, the following discussion focuses on how to adjust the phase state.

4.2. Control of the Phase State

By changing the injection temperature at different injection rates, the relationship between the well bottom temperature and injection time is presented in Figure 6.

Figure 6a–d shows that during the early time of injection, the well bottom temperature changed slightly as the impact of low-temperature CO₂ did not influence the bottom. After an injection period, the bottom temperature decreased considerably. Finally, the temperature became steady. Additionally, with the increase in injection temperature, the steady bottom temperatures shown in Figure 6a–d also increased. Therefore, the high injection temperature could help obtain the supercritical state.

Figure 6a,b shows the temperature distribution when the injection rate was 4 and 2 m³/min, which represents the high displacement rates in CO₂ fracturing wells. It presents that when the injection temperature was −20, −10, and 0 °C, the CO₂ injected into the formation was in the liquid state. When the injection temperature was 10 and 20 °C, the bottom CO₂ was in the supercritical state. If the injection rate was reduced to 0.5 m³/min, as presented in Figure 6c, the bottom temperature increased compared with Figure 6a,b. As shown in Figure 6c, when the injection temperature was −20 °C, the bottom CO₂ was in the liquid state, while the bottom CO₂ was in the supercritical state with the injection temperatures of −10, 0, 10, and 20 °C. If the injection rate was reduced further to 0.1 m³/min, as shown in Figure 6d, the injection temperature had little impact on the bottom temperature, and the bottom temperatures at any injection temperature exceeded the critical temperature. It can be concluded that the low injection rate could increase the wellbore-bottom temperature and help obtain the supercritical state. When the injection rate was relatively small, the injected CO₂ flowed slowly in the wellbore and exchanged heat sufficiently with the surrounding formation. Thus, the supercritical state could be obtained easily by the low injection rate.

Consequently, the supercritical state at the well bottom can be achieved by reducing the injection rate or increasing the injection temperature. To open the reservoir, the injection rate in fracturing wells cannot be reduced considerably. Thus, increasing the injection temperature is an available way to obtain the supercritical CO₂. For CO₂ EOR wells, the injection rate is relatively small. In such cases, the temperature of the wellbore fluid is primarily controlled by the surrounding formation temperature. Therefore, the CO₂ in EOR wells can easily maintain the supercritical state as long as the formation temperature is high enough.

4.3. The Impact of Density Variability on the Flow

As discussed above, the CO₂ injection may reduce the challenge of phase transition between liquid and supercritical states. In the prediction model, the phase transition was coupled by considering the variability of CO₂ physical properties, meaning that the CO₂ density, viscosity, thermal conductivity, capacity, and Joule–Thomson coefficient varied with temperature and pressure. Overall, the variability of density warrants more attention as it denotes the flow compressibility. If the compressibility can be ignored, the model could be simplified. In this section, the influence of density variability caused by the phase transition is discussed.

4.3.1. Criteria of Flow Compressibility

According to the hydrodynamics, the flow compressibility should be determined by the following two conditions. If the two equations are satisfied, the flow can be assumed to be incompressible [29].

$$V\ll c$$

(18)

$$\frac{l}{c}\ll \tau$$

(19)

In Equation (18), c is the sound velocity and V denotes the CO₂ velocity in the tubing. Equation (18) means that fluid velocity should be much less than sound velocity. This relation is satisfied in the injection wells as the fluid velocity in the wellbore is less than the sound velocity. In Equation (19), τ and l are the characteristic time and length, respectively, whose magnitude is the magnitude of time and length when fluid velocity changes significantly. The physical meaning of Equation (19) is that the time when the fluid velocity changes significantly is much greater than the time when sound transits the characteristic length. If the injection rate is relatively small (in CO₂ EOR wells), the variation of fluid velocity could be very small, and the characteristic time is relatively long. Therefore, Equation (19) is satisfied. However, if the injection rate is relatively high (in CO₂ fracturing wells), Equation (19) may not be satisfied as the fluid velocity may change significantly and the characteristic time could be relatively short. Thus, the incompressible assumption may be inapplicable when the injection rate is high. This is the qualitative analysis of compressibility in CO₂ injection wells. The quantitative analysis is presented in the following discussion by error analysis.

4.3.2. Deviations between Incompressible and Compressible Flow

The W-16 well was used to compare the simulated results between incompressible flow and compressible flow. By changing the injection rate, the maximum relative errors between incompressible flow and compressible flow are presented in Figure 7.

Figure 7 illustrates that the relative errors of temperature profiles were small at different injection rates. The maximum error of temperature caused by the incompressible flow was 2.19%. Regarding the pressure distributions, the relative errors were small when the injection rate was less than 2 m³/min. However, when the injection rate was more than 2 m³/min, the relative errors of the pressure profiles were considerable. The maximum error reached 22.32% at the rate of 4 m³/min, which was caused by ignoring the density variability. The high injection rate at the wellhead could lead to high-velocity profiles in the wellbore, so the friction pressure drop was significant. In this case, the pressure distribution was sensitive to the variability of fluid density. Thus, the errors of pressure between compressible and incompressible flow were remarkable when the injection rate was high. Therefore, the compressibility of CO₂ could not be ignored when simulating pressure profiles accurately when the injection rate was relatively high.

In conclusion, to predict the temperature and pressure distribution in low-injection rate wells, the variability of density caused by the phase change can be ignored. However, in fracturing wells whose injection rate is high, the variability of density should be considered to accurately simulate the pressure profiles.

5. Conclusions and Suggestions

With increasing applications of CO₂ in petroleum engineering, the accurate simulation of temperature and pressure in CO₂ injection wells is helpful for effective CO₂ utilization. This work proposed a method to simulate the CO₂ flow in the wellbore. The control of the phase state and the variability of CO₂ density caused by phase transition were analyzed. The main conclusions and suggestions concluded from this work are presented in the following:

• A transient prediction model of CO₂ injection wells was developed, which can simulate the temperature and pressure distributions by the finite difference method. The model was validated using field data.
• The phase state distribution was primarily determined by the wellbore temperature. The phase transition between the liquid and supercritical state may occur during the injection period.
• The supercritical state of CO₂ can be achieved by reducing the injection rate or increasing the injection temperature. For fracturing wells with high injection rates, increasing the injection temperature is possible for the supercritical state. For CO₂ EOR wells with small injection rates, the supercritical state is easily achieved by sufficient heat exchange with the formation.
• When the injection rate is small, the compressibility of CO₂ can be ignored. However, if the injection rate is high, the variability of CO₂ density cannot be neglected as it could lead to significant errors in pressure profiles.

The simulation of CO₂ flow in the wellbore is a challenging task due to its complex conditions. The prediction model could be improved in the following aspects in the future, which was not discussed in detail in this work.

• In fracturing wells, the CO₂ is injected into the formation fractures at the well bottom. Therefore, the influence of the fractures on the CO₂ flow should be studied in the future.
• In some cases, the CO₂ may be injected with water to form CO₂ foam. The behavior of the two-phase flow of CO₂ and water in the wellbore warrants further exploration.