Abstract
CO2 is now considered as a novel heat transmission fluid to extract geothermal energy. It can achieve both the energy exploitation and CO2 geological sequestration. The migration pathway and the process of fluid flow within the reservoirs affect significantly a CO2 plume geothermal (CPG) system. In this study, we built three-dimensional wellbore-reservoir coupled models using geological and geothermal conditions of Qingshankou Formation in Songliao Basin, China. The performance of the CPG system is evaluated in terms of the temperature, CO2 plume distribution, flow rate of production fluid, heat extraction rate, and storage of CO2. For obtaining a deeper understanding of CO2-geothermal system under realistic conditions, heterogeneity of reservoir’s hydrological properties (in terms of permeability and porosity) is taken into account. Due to the fortissimo mobility of CO2, as long as a highly permeable zone exists between the two wells, it is more likely to flow through the highly permeable zone to reach the production well, even though the flow path is longer. The preferential flow shortens circulation time and reduces heat-exchange area, probably leading to early thermal breakthrough, which makes the production fluid temperature decrease rapidly. The analyses of flow dynamics of CO2-water fluid and heat may be useful for future design of a CO2-based geothermal development system.
1. Introduction
The enhanced geothermal system (EGS) is defined as an engineered reservoir that has been created to extract economical amounts of heat from geothermal resources of low permeability and/or porosity [1]. As part of an effort to reduce atmospheric emissions of carbon dioxide (CO2), a novel concept of operating the EGS using CO2 instead of water as the working fluid (CO2-EGS) and achieving simultaneous geologic sequestration of CO2 has been proposed and evaluated [2, 3].
In recent years, a similar concept, the so-called CO2 plume geothermal (CPG) system, has been proposed [4]. The CPG system utilizes existing, naturally porous, high permeability geologic formations (reservoirs) for geothermal energy recovery. The major benefit of the CPG system over the EGS is that the CPG system does not require hydrofracturing, which helps increase fracture permeability but may induce seismicity. The EGS has encountered considerable unfavorable conditions and sociopolitical issues (resistances). Consequently, the CPG system that can use the CO2 sequestration site to recover geothermal energy may be practical.
The advantages of carbon dioxide as a working fluid compared to water in geothermal energy recovery include () large expansivity and compressibility, which result in great density difference between the injector and the producer and reduce the power consumption for circulation on account of buoyancy force, () low viscosity coefficient, which leads to the larger flow rate under fixed differential pressure, () low salt solubility, and () low chemical reactivity; CO2 does not mostly tend to react with rocks. Pruess [3] has built a homogeneous five-spot “fully developed” fractures reservoir model and indicated that, for a given pressure difference between injection and production wells, CO2 would generate 50% larger net heat extraction rate and four times larger mass flow rate compared to water. Atrens et al. [5] drew a conclusion without consideration of frictional pressure that CO2 thermosiphon could produce similar amount of electric power with the same quantity of heat extraction compared to water but with simpler surface equipment. Xu et al. [6] performed wellbore-reservoir coupled five-spot model and pointed out that the specific enthalpy change of CO2 in the wellbore is so small relative to its intrinsic specific enthalpy that the wellbore flow can be regarded as isenthalpic process, and the temperature of the fluids in production well will decrease rapidly with pressure. Compared to water, the advantage of CO2 as a working fluid is more noticeable in low permeability and low temperature reservoirs.
Actually, the investigations mentioned above were carried out under an assumption that the target reservoirs are homogeneous. However, the permeability and porosity heterogeneity of a reservoir significantly affect the CO2 migration and production. Garapati et al. [7] analyzed the effects of multilayered reservoirs on a CPG system; they found that the produced CO2 mass fraction is dominated by high permeability layers and their position within reservoirs. In addition, the heat extraction rate comes down as the permeability of the bottom reservoirs decreases. Yang et al. [8] considered the effects of permeability and porosity heterogeneity on the liquid invasion in tight gas reservoir and drew a conclusion that both the liquid invasion depth and invasion rate increased as the heterogeneity coefficient of permeability decreased. Tian et al. [9] established a two-dimensional model with the consideration of heterogeneity in hydrological parameters of a caprock and indicated that the hydrological heterogeneity significantly affects the containment of intruded CO2 within a caprock. Bu et al. [10] constructed a 2D numerical model to study the influence of high-porosity and high permeability faults within a reservoir and found that the high-porosity and high permeability faults influence CO2 migration and spatial distribution and then increase the CO2 storage capacity of reservoir.
In this paper, we built a three-dimensional geological model with spatially varying permeability and porosity, based on Qingshankou Formation in Songliao Basin, China, to evaluate the effects of hydrological heterogeneity in terms of permeability and porosity on the performance of a CPG system. The temperature change, temporal and spatial distribution of CO2, rate of production flow, heat extraction rate of system, and storage capacity of CO2 were chosen as the metrics to evaluate the performance.
2. Problem Setup
A great deal of detailed information is required to assess the feasibility of injected CO2 as a heat transmission fluid at any specific site and to develop engineering designs for CO2-based geothermal systems. Before moving into site-specific investigations, general features and issues should be explored. This can be done by investigating deep brine systems that abstract site-specific features and thereby attempt to represent characteristics that are common to many such systems. Here, geological characteristics and thermophysical conditions are mainly extracted from the central depression of Songliao Basin, Northeastern China. The basin has a pretty high geothermal gradient and heat flow among sedimentary basins in China and can meet the temperature requirement for geothermal development.
2.1. Conceptual Model and Boundary Condition Settings
To investigate processes of fluid migration and heat exchange within both wellbores and geologic formation, a three-dimensional conceptual model with a size of 10000 m × 10000 m × 100 m, including two wellbores and a reservoir formation (Figure 1), to represent a 100 m thick sandstone layer extending 10 km horizontally in Qingshankou Formation in Songliao Basin, was used. The distance between the injection and production wellbores is 600 m, and both wellbores vertically perforate the whole thickness of reservoir. The length of the model set to be 10000 m laterally is aimed at avoiding the influence of lateral boundaries. Moreover, for the lateral boundaries, Dirichlet-type condition (here, constant pressure and temperature) was imposed. As for the top and bottom boundaries, no-flow conditions were assigned (see the following section for more details). The wellbores with the surrounding formation are also set to be no-flow conditions but the heat exchange is allowed. The mesh of the modeling domain is generated using TOUGHVISUAL [11], which is able to create regular and irregular grids for TOUGH family codes [12]. During simulation, the accuracy and computational demand should be considered. Therefore, the mesh for the area around wellbores is refined, as shown in Figure 1.
2.2. Geological and Thermophysical Parameters
The top of Qingshankou Formation in Songliao Basin is at 2510 m depth and the thickness of this formation is about 250~550 m, including about 100 m thick sandstone layer with porosity of 10%~30% and permeability of 3 to 300 millidarcys which is the main oil and gas reservoir in this formation and with overlying and underlying mudstone as the caprock and bedrock [13]. The average geothermal gradient of the basin is 38.7°C/km [14]. In our model, the depth of the sandstone reservoir is 3000 m with the thickness of 100 m. The initial temperature of the reservoir was set to be 120°C according to the average geothermal gradient and the local mean annual surface temperature (about 4°C). The initial pressure of the model was obtained by hydrostatic equilibrium, and the reservoir was initially saturated with resident water (i.e., CO2 saturation Sg = 0). The details of reservoir geological and thermophysical conditions used in our model are listed in Table 1, and the specification of physical parameters of the wellbores is summarized in Table 2.
2.3. Heterogeneity Implementation
From a realistic view, natural aquifer is intrinsically heterogeneous. The fluids are more likely to flow through the highly permeable zone, especially for CO2, because of its strong mobility. It may lead to an early breakthrough. Among the properties of reservoir, permeability and porosity have the largest influence on flow field and heat extraction. Therefore, we just take them as the main variables to study heterogeneity effects.
As shown in Figure 1, the reservoir extends 10 km laterally, including a one square kilometer subdomain in the center. Here we assume that the subdomain is a geothermal reservoir with heterogeneities in permeability () and porosity (). Neuzil [15] observed that there exists a log-linear relationship between permeability and porosity in argillaceous sediments. Tian et al. [9] drew a similar conclusion by regression analysis. In addition, Nelson [16] suggested a linear relationship between and in quartz sandstone. In this study, we assumed that the above discussed relationship between permeability and porosity is suitable for the Qingshankou Formation, and the relationship can be expressed aswhere and are data-fitted parameters, is the permeability in m2, and is the porosity. Based on the data in Table 1, we figured out the value of and in our model as 10 and −15.52, respectively.
We assumed that the permeability of reservoir is subject to a lognormal distribution, and thereby the porosity obeys a normal distribution inferred from (1). More realistically, the permeability obeys spatial random distribution in the reservoir rather than simply random distribution in probability and statistics. So we introduced the variation function, a concept derived from geostatistics, to depict the spatial distributions of permeability and porosity in the reservoir. In this study, a variance of 0.3 and correlation length of 300 m were chosen [9]. The simulator T2well provides a flexible function inherited from TOUGH2 V2 for permeability modification for each individual grid expressed aswhere is the absolute permeability of grid , as specified in data block ROCKS and is the permeability modification coefficient which can be defined internally or externally and may be provided as part of the geometry data in block ELEME and then employed to multiply the absolute permeability for each grid in the subdomain. When the permeability modification is accomplished, the porosity is yielded according to (1) and stored in block INCON. More details about permeability modification coefficient generation and heterogeneity realization can be found in previous research [9]. The permeability and porosity range is shown in Table 3.
Infinite types of heterogeneous permeability and porosity fields can be generated randomly. As shown in Figure 2, six cases with heterogeneity in permeability and porosity were selected to represent the conditions: (i) high permeability belt connecting two wellbores (hh1, hh2); (ii) highly permeable zones expanding away from the production wellbore (hp1, hp2); and (iii) low permeability zone between the two wellbores (hl1, hl2). They are not expected to match with any existing geologic formation; we aim to make some theoretical analysis and quantitatively evaluation of the effects of heterogeneities on a CPG system and elucidate the affecting mechanism.
(a) hh1
(b) hh2
(c) hp1
(d) hp2
(e) hl1
(f) hl2
2.4. Modeling Scenario
The CO2 migrates simultaneously as the CO2 is injected into the reservoir through injection wellbore, which pushes the existing resident water to move. Therefore, the production fluid is expected to be water at the early age. So the injection and production pressures at the wellhead are set to be 12 MPa and 0.1 MPa (atmospheric condition) to make CO2 get breakthrough as soon as possible. After a period of propagation, CO2 passes through the reservoir and can be observed from the production wellhead, and then the saturation of CO2 (Sg) in production fluids increases rapidly, and the temperature of carbon dioxide decreases sharply with the pressure decline in the production wellbore, which is called the Joule-Thomson effect. To get a higher production fluid temperature and control the flow rate for a stable operation of the binary system, a constant pressure of 8.0 MPa is imposed both on the wellhead of production and injection after the production well is occupied by CO2 (about 3 years). CO2 is injected into the reservoir at a constant temperature of 15°C for 20 years.
3. Simulation Approach
3.1. Governing Equations
The process of fluid flow in subsurface can be divided into two parts, flow in the wellbores and flow in the reservoir formation. Therefore, we employed the wellbore-reservoir coupled simulator T2well [17–19]. It is an extension of the general multiphase, multicomponent, nonisothermal simulator TOUGH2 V2 [12]. The program assigns the wellbore and reservoir as two subdomains, in which flows are controlled by appropriate laws, respectively. In the reservoir, the flow is described as Darcy’s law and in wellbore, it is momentum conservation. To consider a comprehensive description of the thermophysical properties of H2O-CO2 mixtures, the ECO2N V2.0 module was introduced, which reproduces fluid properties including density, viscosity, and specific enthalpy, largely within experimental errors under the temperature and pressure conditions of 10°C < < 300°C, < 600 bar [20]. These fundamental flow equations, used in the T2Well code [18], are summarized in Table 4.
In Table 4, means the accumulation term, and represents the flux term of mass or energy. represents phase and is the index for fluid species. For mass conservation, are density, volumetric fraction, and velocity of phase-, respectively. is the mass fraction of component in phase-. For energy conservation, , , and are heat conductivity, specific enthalpy of fluid, and inclination angle of wellbore, respectively. For momentum equation, is perimeter of wellbore, and is the wall shear stress.
3.2. Drift-Flux Model
Drift-Flux Model [21, 22] is introduced to calculate the two-phase velocities of CO2-water mixtures in wellbores by the following equations.
First, the velocity of gas phase can be described by the constitutive relation as below: where the profile parameter is used to account for the effect of local gas saturation and velocity profiles over the pipe (wellbore) cross-section [23]; is the volumetric flux of total mixture, which is defined as Thus, the velocity of liquid could be determined asThen the mixture velocity can be calculated in light of By inserting (3) and (5) into (6), the volumetric flux could be described as a function of mixture velocity and the drift velocity , as (7) shows: where , and it is the profile-adjusted average density. The major task is now to calculate the mixture velocity and the drift velocity.
4. Results and Discussion
4.1. Migration of CO2 in the Reservoir
The distribution of CO2 saturation in the reservoir after 3 years of injection is shown in Figure 3. By comparing Figure 3(a) through 3(f) (6 heterogeneous cases) with the homogeneous cases of Figure 3(g), it can be found that the migration of CO2 is significantly affected by the media (permeability and porosity) heterogeneity. In the homogenous reservoir, the CO2 saturation distributes circularly around the injection well. However, in the heterogeneous cases, the CO2 tends to flow through the highly permeable zone and leads to a preferential flow. It bypasses the low permeability zone to reach the production well, even if the flow path is longer (shown in Figures 3(e) and 3(f)). The positions holding high CO2 saturation match well with highly permeable zone. For the cases that the highly permeable zones locate far away from the production well (Figures 3(c) and 3(d)), there would be more injected CO2 stay in reservoir instead of being extracted from the production well.
(a) hh1
(b) hh2
(c) hp1
(d) hp2
(e) hl1
(f) hl2
(g) Homo
4.2. CO2 Saturation in Production Fluids
Being injected through the injection wellbore, CO2 displaced the original water in the reservoir formation, and then water was first produced from the production well. At early stage, the production fluid was pure water, and CO2 plume expanded within reservoir over time. When it reached the production well, the output fluids became into a mixture of water and CO2. As can be seen in Figure 4(b), the time of CO2 breakthrough was quite different in different cases. In homogeneous case, it takes about 1.3 years for CO2 to reach the production well, while the breakthrough time of most heterogeneous cases is earlier; for example, it only takes 0.8 years for CO2 to get through the reservoir and reach the production well in hh2 case, within which there is a highly permeability zone between two wells. Correspondingly, in hl1 and hl2 cases, the breakthrough time is 1.4 years and 1.5 years, respectively, due to the low permeable belt between the injection and production wells. After CO2 reached the production well, the saturation of CO2 in the output fluids increased rapidly to over 0.9 (Figure 4(a)). The dive of curves in the figures is caused by the increasing of production pressure (from 0.1 MPa to 8 MPa as mentioned before). After that, the saturation of CO2 recovered smoothly and stayed at about 0.92 during the rest operation time.
(a)
(b)
4.3. Injection and Production Rate of Fluids
The evolution of injection rate of CO2 of different cases is presented in Figure 5. As shown in the figure, the injection rate of CO2 rises rapidly after it reaches the production well. When injection and production wellhead pressures are imposed to be 8 MPa, the cyclic pressure difference () reduces to zero from 11.9 MPa, the injection rate of CO2 instantly fell and stayed relatively constant during the rest time of simulation. The injection rate of CO2 in hh1 and hh2 is relatively higher compared to other cases due to its high permeability and porosity zone between the two wellbores.
Figure 6(b) shows the evolution of production rate of water. It can be seen that the rate increases dramatically like an eruption, when CO2 appears from the production wellhead. This phenomenon is caused by the CO2 accumulation in production; when it exceeds a certain quantity, an eruption event would take place. This phenomenon can likely account for some cold geysers eruptions [24, 25]. After the eruption, the content of water in the mixture fluids goes down. At the instant of raising the production pressure and reducing the injection pressure to 8 MPa, both the production rate of water and CO2 fall off. Note that the production rate of water declines continually but that of CO2 rises gradually. Even though the heat transmission fluid is CO2, water still takes a large proportion of the production fluids because of the huge density difference in production well. Simulation results show that low permeability zone between the wellbores (hl1 and hl2) does not correspond to the lowest production rate of fluids as we expected. On the contrary, water production rate is relatively larger compared to other cases. This is because the resident water saturation at production well bottom in hl1 and hl2 is higher than that of other cases. On the other hand, the highly permeable belts connecting two wells may cause the high production rate of CO2 and low production rate of water (hh1 and hh2). As can be seen in the graph, the fluids production rates of hp1 and hp2 are really low during the simulation time, which could affect the economic feasibility of the CPG system.
(a)
(b)
4.4. Temperature of Production Fluids
The temperature of production fluids is affected by the flow path and heat exchange within the reservoir. As shown in Figure 7, the temperature of production water distributes between 100 and 115°C before CO2 reaches the production wellhead, which is somewhat less than the reservoir initial temperature (120°C). This temperature drop in the production well is due to the heat exchange (loss) with the surrounding rock around the wellbore. Consequently, when CO2 begins to produce in the wellhead, the temperature of mixture output fluids drops rapidly. It is caused by the Joule-Thomson expansion in the wellbore with the pressure decrease from the downhole to the wellhead. It can be seen that the temperature drop of output fluids from the downhole to wellhead is approximate 40°C. Enhancing the production pressure increases the temperature of output fluids by 10 to 14°C.
Comparing the temperature of production fluids between different cases, it can be found that, to a certain extent, the temperature of output mixture fluids is proportional to the water output rate. It can be explained from two aspects. Firstly, water can counteract the Joule-Thompson effect and help keep the temperature. Secondly, higher water output rate implicates that more space of reservoir is occupied by CO2, and more energy is extracted by the output fluid. That is just the reason why the cases with low permeability zone between the two wells correspond to a higher output fluids temperature (hl1 and hl2).
4.5. CO2 Storage
To investigate the CO2 storage capacity in the reservoir, we define the storage rate of CO2 as the injection rate minus production rate of CO2 (). The storage rate of CO2 of different cases is shown in Figure 8(a). More CO2 is stored in the cases with highly permeable zone located away from production well. The storage amount of CO2, that is, the storage rate of CO2 integral of time, is shown in Figure 8(b). Similarly, the storage amount of CO2 of cases hp1 and hp2 (with highly permeable zone deviating from production well) is relatively higher compared to other cases. However, the storage amount of CO2 is also affected by cumulative injection of CO2. So the storage ratio of CO2 is calculated by storage amount of CO2 divided by cumulative injection of CO2 and is shown in Figure 9. In cases hh1 and hh2, CO2 tends to migrate to production wellbore through high permeability zone between the two wellbores and the storage ratio of CO2 in the reservoir is less than homogeneous case. On the contrary, in cases hp1 and hp2 (with high permeability belt deviating from producer) a portion of CO2 will transport away from production wellbore and trapped in the reservoir. The storage ratio of cases hp1 and hp2 is 52.7% and 55.3%, respectively, which is relatively higher compared to homogeneous case (45%).
(a)
(b)
4.6. Temperature Distribution of Reservoir
The temperature drop in reservoir is caused by heat transfer between rock and fluids. After 20 years of production, the temperature of the rock matrix between the injector and producer is significantly reduced. It is highly affected by the hydrological heterogeneity. In cases hl1 and hl2 (Figures 10(e) and 10(f)), the temperature drop around the production wellbore is less than other cases due to the low permeability and porosity zone between the two wellbores.
(a) hh1
(b) hh2
(c) hp1
(d) hp2
(e) hl1
(f) hl2
(g) Homo
4.7. The System Heat Extraction Rate
The net heat extraction rate is calculated by the following [26]:where is the net heat extraction rate, is the production flow rate, is the specific enthalpy of injection fluids, and is the specific enthalpy of production fluids.
The simulated heat extraction results (Figure 11) indicate that media heterogeneity affects the heat extraction rate greatly. In cases with low permeability zone between production and injection wells (hl1 and hl2), the heat extraction rates are similar with the homogeneity case and stay relatively constant in the middle and later periods. Cases hp1 and hp2, with high permeability zones extending far away from production well, have low heat extraction rates during the entire simulation time. The other heterogeneous cases with high permeability zone between two wells (hh1 and hh2) obtain the maximum extraction rates at early time (2 to 3 years) and decrease quickly in the later periods.
5. Concluding Remarks
We have built a three-dimensional wellbore-reservoir coupled model with consideration of permeability and porosity heterogeneity based on the geological and thermal-physical conditions of Songliao Basin, China. A total of 7 case simulations were performed. The following conclusions can be drawn.
Heterogeneity of reservoir’s hydrological properties (in terms of permeability and porosity) affects the migration of CO2 in the reservoir significantly. Due to the strong mobility of CO2, as long as a highly permeable zone exists between the two wells, it is more likely to flow through the highly permeable zone to reach the production well, even though the flow path is longer. The preferential flow shortens circulation time and reduces heat-exchange area, probably leading to early thermal breakthrough, which makes the production fluid temperature decrease rapidly. The highly permeable zone located away from production wellbore is more in favor of the storage of CO2 in the reservoir; however, it leads to quick decline in temperature of production fluids and heat extraction rate, which could affect the economic feasibility of the CPG system. Both of cases mentioned above should be paid more attention in the future design of an actual CPG demonstration project.
The range of problems concerning the process in CO2-based geothermal systems is very broad. The present simulation results are specific to the conditions and parameters considered. The “numerical experiments” give a detailed understanding of the dynamic evolution and provide useful insight into fluid flow and thermal dynamic processes along the wellbores and reservoir.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This work was jointly supported by National Program on Key Research and Development Project (no. 2016YFB0600804), by the National Natural Science Foundation of China (Grant no. 41572215), by the 111 project (no. B16020), and by the project of CO2 geological storage in China Junggar basin of China Geological Survey (Grant no. 121201012000150010).