CFD Convective Flow Simulation of the Varying Properties of CO2-H2O Mixtures in Geothermal Systems

Numerical simulation of a geothermal reservoir, modelled as a bottom-heated square box, filled with water-CO2 mixture is presented in this work. Furthermore, results for two limiting cases of a reservoir filled with either pure water or CO2 are presented. Effects of different parameters including CO2 concentration as well as reservoir pressure and temperature on the overall performance of the system are investigated. It has been noted that, with a fixed reservoir pressure and temperature, any increase in CO2 concentration leads to better performance, that is, stronger convection and higher heat transfer rates. With a fixed CO2 concentration, however, the reservoir pressure and temperature can significantly affect the overall heat transfer and flow rate from the reservoir. Details of such variations are documented and discussed in the present paper.


Introduction
Simultaneous power generation and geosequestration make CO 2 a very attractive choice for geothermal power plants. As such, carbon-dioxide-based engineered geothermal systems (CO 2 -EGS) have been previously proposed as an alternative to water-based EGS systems [1]. Subsequent studies added further details and reported possibility of improved energy extraction [2][3][4][5][6][7]. Interestingly, the "dry-out period, " or transition from an initially water-filled EGS system to a CO 2 -rich one [8][9][10][11] has received a lot of attention mainly because of possible mineral dissolution and precipitation as a result of changes in the reservoir fluid composition and reservoir permeability alteration [12][13][14]. Additionally, [15] investigated the effects of CO 2 -rich phase compositions on the production flow rate and the heat extraction from the reservoir. What is yet to be reported in the literature is a detailed numerical simulation of a water-CO 2 mixture filling a reservoir. Simple thermodynamic analysis of a reservoir shows that more heat can be extracted (compared to a watersaturated reservoir) mainly because a CO 2 -water mixture is more buoyant than pure water. This could significantly affect the energy extraction from a water-saturated reservoir. Most of these reservoirs suffer from permeability drops with depths. As such, having a more buoyant fluid which could move up against gravity, at least partially, not only is desirable but also is going to significantly affect the drilling cost which is expected to grow exponentially with the well depth [16][17][18].
The role of fluid migration in the build-up of heat in underground geothermal systems is not well understood. It is known that degassing of CO 2 with isotopic composition indicating mantle-sources in regions of tectonic activity is associated with locally elevated geothermal temperatures [19][20][21]. Convective fluid plumes may play a role in enhancing heat flows from the mantle to geothermal reservoirs and within the reservoirs themselves. Additionally, convection within a geothermal reservoir may enhance the productive life-time of geothermal reservoir systems by enhancing heat supply from underlying strata and by ensuring a more even distribution of thermal energy throughout the reservoir, that is, by off-setting localized cooling along major flow paths. This could also imply that for a water-saturated reservoir which is classified as nonproductive, mainly due to lack of convection currents, injection of CO 2 can lead to formation of convection cells and thus facilitate heat extraction. In what follows, a numerical analysis of this problem is presented to quantitatively measure the improved convective flow patterns and enhanced heat transfer from the reservoir. This is intended to provide insight as to the possible mechanisms by which CO 2 presence in or addition to underground reservoirs could lead to enhancement of convective heat transfer.

Modelling
The reservoir is modelled as a bottom-heated square box with adiabatic lateral boundaries and a cold top wall as Figure 1 shows. The cold and hot temperatures are varied from 331 to 431 K and 416 to 516 K, respectively, in a way that the hot-cold temperature difference remains at 85 K for each case. For constant properties, and of course with the same reservoir size, porosity, and permeability, one would expect the results to be the same as long as the temperature difference is not altered. It will, however, be shown in the forthcoming sections that this is not the case in our problem as properties significantly vary with both temperature and pressure. The reservoir porosity is fixed at 0.05 and the permeability-length product is kept constant at 10 −11 m 3 similar to Haghshenas Fard et al. [7] with no through-flow. The reservoir pressure is varied from 20 to 60 MPa (equivalent to the hydrostatic pressure of approximately 2 to 6 km of water) to cover a wide range of practical applications for geothermal development.
The real (i.e., nonideal) properties are determined for pure CO 2 using a Helmholtz-free-energy based equation of state [22] and for pure H 2 O using the International Association for the Properties of Water and Steam equation of state [23]. The use of a single-phase basis for modelling is to enable examination of the effect of changing properties without additional complexities of two-phase flow separation behaviour. Further comment is provided in the discussion. The reservoir properties are then obtained similarly using the porosity and solidity as the weight-functions as follows: One notes that in the above formulation changes in the reservoir porosity and permeability are not taken into account while fluid properties are updated in each iteration.
An in-house code, used by Hooman and Gurgenci [24], was cross-validated with the commercially available software ANSYS-FLUENT and used to create and mesh the geometry and finally solve the governing equations. The governing equations were derived as for standard two-dimensional heat transfer conditions, based on the assumptions of adiabatic lateral walls and constant temperature for the top and bottom boundaries, and are as follows: subject to the boundary conditions illustrated in Figure 1.

Numerical Details
Grid independence was verified by running the software on different combination of grid sizes. It was observed that the results changed less than 2% when a 100 × 100 mesh system is used instead of a finer mesh with 200 × 200 grid points. Results are also verified for constant property free convection of water in a porous cavity, that is, the Darcy-Benard problem. It was noted that the correlation between Nusselt number (Nu) and Rayleigh number (Ra) Nu = Ra/40 best fits out numerical data, as Figure 2 shows. Ra in this instance is determined as per As a further check on the accuracy of our results, variable property Darcy-Benard free convection of pure water in a porous cavity was investigated to observe that using the reference temperature approach of Hooman and Gurgenci [24] the above correlation can still be used within 5%.

Results and Discussion
In what follows we focus on free convection heat and fluid flow of a water-CO 2 mixture in a porous cavity. We use Nu and maximum flow rate as our metrics to evaluate the The Scientific World Journal strength of convective flow patterns. Nu is the total heat transfer divided by that of pure conduction through the same cavity (no convective flow patterns). As such, any Nu value in excess of unity shows some degree of convection. Obviously, higher Nu values mark stronger convective cells. The flow rate reported here is the one induced by free convection only, that is, without a well-head pump or any other suction/injection mechanisms. We systematically change the CO 2 mass fraction from zero (pure water) to unity (pure CO 2 ) over a range of reservoir pressure and temperature in a way that the hot-cold temperature difference remains the same. For a constant property subcritical fluid flow, one would expect that, with the same temperature difference and, hence, the same Rayleigh number (Ra), the overall heat transfer and fluid flow will not alter. However, as CO 2 is supercritical within the range of conditions of underground reservoir systems, that is not the case for mixtures of CO 2 and H 2 O, as demonstrated by Figure 3. This figure shows Nu versus CO 2 mass fraction at 20 MPa with the same hot-cold temperature difference but with different hot and cold temperatures as denoted on the plots. As seen, the heat transfer increases with CO 2 mass fraction for any given ℎ and combination. Furthermore, moving from pure water to pure CO 2 , the increase in heat transfer is significant; about one order of magnitude is the minimal heat transfer augmentation. More interestingly, however, is the fact that Nu is the highest with the lowest (and obviously lowest ℎ to maintain the same Δ of 85 K) mainly because the lower temperatures are closer to those of pseudocritical conditions where Ra is expected to reach a maximum value; see also Forooghi et al. [25][26][27][28]. This is obviously in favour of low temperature geothermal  Figure 3. That is, heat transfer increases for temperatures close to pseudocritical conditions. Moreover, based on plots in the same chart, increasing the pressure leads to lower heat trasnfer rates for a fixed CO 2 mass fraction and temperature. This could be explained as the obvious decrease in compressibility and increase in the fluid density with higher pressures, with a fixed fluid temperature, which will lead to lower thermal expansion coefficients. As a result, at the same temperature, either of the two fluids will be less buoyant at higher pressures when compared to lower ones, so will be the mixture in the absence of any phase transitions.
The convective flow rates also reflect a dependence on compressibility, as demonstrated by Figure 5. The dimensionless flow rate (normalized stream function on the vertical axis) is obtained by normalizing the actual flow rate with appropriate scales for velocity, area, and density: Mathematically, it means that we used the group to nondimensionalize the flow rate. It needs to be mentioned that the choice of these parameters is optional but we tried to use constant values for density and length to make it easy for the reader to generate estimates, based on our calculation, for expected flow rates through a given reservoir. Moreover, what we are more interested in is the trend of the flow rate plot against the mass fraction than the actual flow rate values. In doing so, the (constant) density of water at atmospheric condition is used where the unit area is used defined as the length of the cavity multiplied by unity (1 m). The flow velocity, for single-phase constant property case, is assumed to be linearly proportional to the product of the thermal diffusivity and Ra 1/2 and inverse linearly proportional to the cavity length; for example, The flow rate is given bẏ The Scientific World Journal The product of thermal diffusivity and density is independent of density and leaves us with a group / . Consequently the mass flow rate scale is represented aṡ While and are calculated under standard atmospheric conditions, Ra is affected by fluid property variation and following the use of (1)- (2). This flow rate here is the buoyancy-induced flow rate due to changes in fluid density. The Scientific World Journal normalized by using (8) above are presented in Figure 5 for different mass fractions, pressure, and temperature combinations. Similar to Nu plots, one notes that the mass flow rate is sensitive not only to the temperature difference but also to the actual wall temperature values. Furthermore, higher CO 2 mass fraction leads to higher flow rates. It can be noted, moving from Figures 5(a) to 5(c), that flow rates are less sensitive to pressure as ℎ is increased. With a fixed ℎ and , one notes different trends in flow rate when pressure changes. Depending on the temperature values, an increase in pressure can either increase ( Figure 5(a)) or decrease the flow rate ( Figure 5(b)). This work is an initial analysis of the role of CO 2 enhancement of convective heat transfer within geothermal reservoirs. It deliberately assesses the behaviour of a single-phase mixture of the two components. Further work is necessary to extend this to account for multiple phases. There are three particular qualitative effects through which multiphase flow is expected to alter the results presented here.
(1) Transient exsolution of dissolved CO 2 as bubbles within the two-phase region should lead to local enhancement of convective flow around the bubble due to its upwards buoyancy-driven motion. One expectation of this would be an increase in the gradient of convective heat flux with mass fraction (i.e., / ) at the bubble line, where CO 2 begins to exsolve from the H 2 O phase.
(2) Relative permeability within the two-phase region would act to reduce enhancement of flows, as the reservoir permeability to the minor phase within a two-phase flow is typically substantially reduced.
These expected qualitative behaviours require further analysis accounting for multiphase flow behaviour to determine their relative contribution to overall convective heat flow enhancement. However, the sum of these changes is not expected to reverse the overall trend demonstrated here, of increased convective heat flux as CO 2 is added to the reservoir system. Considering that the results presented here indicate that CO 2 may enhance flow rates by up to a factor of 2.67, we conclude that this is a potentially important mode of heat transport within geothermal reservoirs and warrants further study. We anticipate the next steps to be consideration of the additional flow behaviours when multiple phases are present.

Conclusion
The effect on convective heat transport within a closed reservoir system of varying fluid composition was analysed by CFD modelling of a single-phase fluid with properties derived from a composition-dependent average of pure CO 2 and pure H 2 O. As compositional properties were varied from H 2 O toward that of CO 2 , substantial increases were observed in Nusselt number (by a factor of up to 10) and normalised stream function (by a factor of up to 2.67). We conclude that this indicates substantial increase in convective heat transport.
Convective heat transport may be further modified by multiphase flow behaviours, and we conclude that, due to the potential magnitude of heat flow enhancement by the addition of CO 2 , further research exploring the effect of these behaviours on heat transport is warranted.
We anticipate this finding to have implications for the study of natural geothermal reservoirs, where the role of dissolved gas exsolution on heat transfer enhancement remains unquantified. Additionally, these findings may be of potential interest with regard to CO 2 injection into geothermal reservoirs, as it may lead to improved productivity through the mechanisms elucidated in this work.
The magnitude of this effect on the thermal productivity of geothermal power plants is difficult at this stage to quantify and probably not meaningful to speculate on due to the limitation of assessing only single-phase flow behaviours within the reservoir. These findings do demonstrate, however, that reservoirs with elevated CO 2 content will experience greater convective heat transfer and therefore be of comparatively higher temperature (and therefore of greater resource value). Additionally, it can be concluded that any increase in CO 2 content of an existing reservoir will enhance convective flow behaviours (although the true magnitude of this effect will depend on two-phase flow behaviour as well as the particulars of the reservoir) and consequently will enhance the productivity and/or the longevity of geothermal energy extraction.