- Email: [email protected]
Numerical analysis of heat extraction performance of a deep coaxial borehole heat exchanger geothermal system
Accepted Manuscript Numerical Analysis of Heat Extraction Performance of a Deep Coaxial Borehole Heat Exchanger Geothermal System
Xianzhi Song, Gaosheng Wang, Yu Shi, Gensheng Li, Ruixia Li, Zhengming Xu, Rui Zheng, Yu Wang PII:
S0360-5442(18)31588-3
DOI:
10.1016/j.energy.2018.08.056
Reference:
EGY 13532
To appear in:
Energy
Received Date:
09 April 2018
Accepted Date:
07 August 2018
Please cite this article as: Xianzhi Song, Gaosheng Wang, Yu Shi, Gensheng Li, Ruixia Li, Zhengming Xu, Rui Zheng, Yu Wang, Numerical Analysis of Heat Extraction Performance of a Deep Coaxial Borehole Heat Exchanger Geothermal System, Energy (2018), doi: 10.1016/j.energy. 2018.08.056
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT Numerical Analysis of Heat Extraction Performance of a Deep Coaxial Borehole Heat Exchanger Geothermal System Xianzhi Song1*, Gaosheng Wang1, Yu Shi1, Gensheng Li1, Ruixia Li2, Zhengming Xu1, Rui Zheng1, Yu Wang1 (1. State Key Laboratory of Petroleum Resources and Prospecting, China University of Petroleum, Beijing, Beijing 102249, China; 2. Sinopec Star Petroleum Co. Ltd, Beijing 100083, China)
ABSTRACT: The coaxial borehole heat exchanger (CBHE) can effectively develop deep geothermal resources with larger heat exchange area than traditional U-tube BHE. However, few studies focus on its parameters optimization and temperature field dynamic response. Hence, an unsteady-state heat transfer model is established for a deep CBHE. The finite difference method is adopted to solve the model, and the model is validated with experiment data. The effects of key factors on the CBHE performance are comprehensively analyzed. Simulation results depict that the outlet temperature decreases greatly at the initial stage, and then remains relatively stable. There is a critical value for flow rate to obtain higher thermal power with appropriate pressure drop. The insulated pipe can significantly reduce heat dissipation of working fluid in central pipe. There is an optimal insulation length in terms of heat insulation effect and costs. The thermal conductivity of cement sheath has remarkable influence on the thermal process and partial cementing with high thermal conductivity cement is favorable. The temperature field analysis shows that the impact scope is limited, and the temperature in the reservoir before the next heating period would almost recover its original condition, verifying the sustainability of CBHE system. KEY WORDS: geothermal energy; coaxial borehole heat exchanger; unsteady-state heat transfer; numerical simulation, heat extraction performance 1. Introduction The massive consume of the fossil fuels have caused many problems, such as resource depletion, global warming and environmental pollution [1, 2]. Therefore, the renewable energy such as geothermal energy is playing an increasingly important role in the district energy supply and environmental improvement [3]. The borehole heat exchanger has been widely used to extract geothermal energy for space heating in the residential and commercial buildings, and the U-tube borehole heat exchanger is the most common form [4, 5]. However, the coaxial borehole heat exchanger (CBHE) is preferred to develop the deep geothermal resources for larger heat exchange area and better performance [6, 7]. Investigation on the heat transfer process of CBHE system has been carried out by many researchers. Roland N. Horne [8] first established a quasi-steady heat transfer analytical model, and studied the effects of flow direction, inlet temperature and outer tube radius. To optimize the flow rate and geometric size, P.J. Yekoladio et al. [9] developed a steady state analytical model, while the dynamic change of temperature field in formation was ignored. Richard A. Beier et al. [10] proposed a transient heat transfer model solved by Laplace transform to investigate the borehole thermal resistance. David Gordon et al. [11, 12] presented composite cylindrical source (CCS) model to analyze the short-term performance of CBHE system, then the influence of pipe materials was evaluated by a semi-analytical model. Analytical models contribute to explore the rigorous relationship of each variable, while the numerical solution is also a convenient tool to study the performance of CBHE. Morita et al. [13, 14] employed explicit finite-difference method to solve the heat transfer model, and investigated the thermal characteristics. Henrik Holmberg et al. [6] discretized the energy equation by implicit finite difference scheme, and results showed that the * Corresponding author. E-mail address: [email protected] (X. Song).
ACCEPTED MANUSCRIPT deeper boreholes was beneficial for heat extraction. Then Mohamed Salah Saadi et al. [15] built a transient 2-D FVM (Finite volume method) model considering the seasonal effects for CBHE. Previous studies have made significant contributions to comprehend the thermal process of CBHE. However, few studies focus on the parameters optimization, for example the critical values of flow rate, thermal insulation length, cementing length, etc. And the investigation on the temperature field dynamic response in reservoir is still lacking, especially after a recovery period, which is vital for well spacing determination and sustainability analysis. Additionally, the commercial simulation software, such as COMSOL and FEFLOW, also provides a new platform to study the thermal process. Zanchini et al. [16, 17] established a 2-D model through COMSOL, then found that the thermal conductivity of cement played a remarkable role in the heat transfer. Xianzhi Song et al. [18, 19] presented a 3-D transient state model coupled with groundwater flow by COMSOL, and evaluated the heat transfer characteristics of different borehole heat exchangers. Based on FEFLOW code, some researchers analyzed the subsurface thermal environmental changes and seasonal heat storage efficiency [20, 21]. Nevertheless, the mesh generation for many thin layers of large scale model, such as casing and cement sheath, is not convenient in commercial software. Besides, Guodong Cui et al. [22] proposed a promising CBHE with single horizontal well, but the vertical CBHE is recommended in view of drilling cost and complexity. To evaluate the actual heat extraction performance of CBHE, a field experiment was carried out at HGP-A well in Hawaii, which also provided reliable data for model validation [14, 23]. In 1993, a deep CBHE in Weissbad, Switzerland, was operated commercially, but its outlet temperature was lower than expected. T. Kohl1 et al. [24, 25] attributed this to the inner tube without insulation and the poor contact between casing and formation. Furthermore, a distributed thermal response test (DTRT) for CBHE was conducted to accurately analyze borehole thermal resistances [26]. China has the most installed thermal capacity for geothermal direct-use, but its geothermal power development is far below than that of the top 10 countries, such as United States, Mexico and Iceland [27-29]. The main reason is the relatively low temperature in single geothermal well, which leads to a limited heat output. Thus, there is an urgent need to exploit the deep geothermal resources in China. And the CBHE with a vertical wellbore would be employed to extract geothermal energy in Xiong'an New Area for its better performance and cleanliness. However, there are few comprehensive researches of parameters optimization, including flow rate, inlet temperature, cementing length, insulation section length, especially for a deep CBHE. And it is still lack of further analysis on the dynamic change of temperature field in formation. Hence, an unsteady-state heat transfer model is proposed in this paper. Then the finite difference method is adopted to solve the mathematical model, and the model is validated with experiment data. The main findings can offer guidance for the optimal design of CBHE system. 2. Model establishment 2.1 Physical Model The CBHE system consists of high-pressure pump, vertical wellbore and geothermal heat pump in Fig. 1. A working fluid is injected from the annulus, and then produced through central pipe. Forced convection occurs in the wellbore, while the thermal process in solid simply involves heat conduction. Through heat pump, the geothermal energy carried by the high-temperature fluid can be used for space heating or power generation. Additionally, because of the high-pressure environment in a deep well, the conventional BHE materials with low thermal conductivity [6, 30],
ACCEPTED MANUSCRIPT such as polyethylene and polypropylene, are difficult to meet the strength requirement. Hence, a three-layer vacuum tube is employed, and the pressure of insulation layer remains 0.2 bar by a vacuum pump, which avoids air circulation.
Fig. 1. The schematic diagram of the flow and heat transfer in CBHE 2.2 Mathematic model Due to the symmetry of simulation domain and boundary conditions, a 2D unsteady-state heat transfer model is set up based on the energy balance relationship among the heat exchange medium. And the circulating fluid, inner tube, casing and cement sheath are simplified into one-dimension, while we adopt two-dimensional grid for the reservoir simulation. And the assumptions in this model include: (1) The heat generation by viscous friction of working fluid is negligible compared to the heat exchange between working fluid and its surrounding medium; (2) Regarding water as the working fluid, and its heat transfer form is limited to forced convection; (3) Only conductive heat transport is considered in formation. 2.2.1 Heat transfer model in central tube According to Fig. 1, the energy balance equation is established for the working fluid. The second term on the left side refers to heat transfer from annular fluid. The thermal process involves forced convection and heat conduction through central tube. The expression is as follows:
m qC pm
T f 1 z
Tf 2 Tf 1 R
mC pm r12
T f 1 t
（1）
where Tf1 and Tf2 (°C) are the temperature of working fluid in central pipe and annulus, respectively. R ((m·°C)/W) is the thermal resistance, which is similar as electric resistance. This variable is defined as:
R (W/(m2·°C))
ln r2 / r1 ln r3 / r2 ln r4 / r3 1 1 2 r1h1 2w1 2w 2 2w3 2 r4 h2
（2）
where h is the convective heat transfer coefficient, which is affected by temperature due to the change of physical properties of working fluid. This coefficient is defined as:
ACCEPTED MANUSCRIPT
h
Nu f D
（3）
where Nu is Nusselt number, which characterizes the intensity of convective heat transfer. It is calculated by Gnielinski equation [31]:
f Re 1000 Pr 8 Nu 1 12.7 f
8
Pr
23
1
（4） where Re is Reynolds number. Pr is Prandtl number. f is Darcy friction factor of turbulent flow in a tube, which is calculated by Filonenko equation [32]:
f
1
1.82 lg Re 1.64
2
（5）
Re=3000~6×106,
The application range of Eq. (5) is: Pr=0.5~2000, L/D>10, which represents the ratio of casing length to diameter. Then the hydraulic diameter is introduced for the flow in annulus:
Deq 2r5 2r4
（6）
where r4 and r5 (m) are the inner and outer radius of the annular space, respectively. 2.2.2 Heat transfer model in annular space The energy balance equation for annular fluid is established based on the thermal process. Compared with Eq. (1), there is one more item on the left side in Eq. (7), which represents the heat transfer from casing, as shown below:
m qC pm
T f 2 z
Tf 1 Tf 2 R
2 r5 h3 Tw T f 2 mC pm r52 r42
T f 2 t
（7）
2.2.3 Heat transfer model in casing Due to the temperature difference, the heat would transfer from cement sheath to casing, which is represented by the third item in Eq. (8). Also, it needs to consider the forced convection heat transfer. The equation is described by:
2Tw 2wc w r r 2 2 r5 h3 T f 2 Tw Tc Tw z ln((r7 r6 ) / (r6 r5 )) 2 6
2 5
wCw r62 r52
Tw t
（8）
where λwc (W/(m·°C)) is the equivalent thermal conductivity of casing and cement, which is determined by:
wc =
ln (r7 r6 ) (r6 r5 ) ln 2r6 (r6 r5 ) ln (r7 r6 ) 2r6
w
c
（9）
where λw and λc (W/(m·°C)) are the thermal conductivity of casing and cement sheath, respectively. 2.2.4 Heat transfer model in cement sheath The heat transfer process in cement sheath is only limited to heat conduction, and its energy balance
ACCEPTED MANUSCRIPT equation is set up as:
2Tc 2wc 2cs c r r 2 Tw Tc Ts Tc z ln((r7 r6 ) / (r6 r5 )) ln((r8 r7 ) / (r7 r6 )) 2 7
2 6
c Cc r72 r62
Tc t
（10）
where λcs (W/(m·°C)) is the equivalent thermal conductivity of cement and formation, and it is identified by:
cs =
ln (r8 r7 ) (r7 r6 ) ln 2r7 (r7 r6 ) ln (r8 r7 ) 2r7
c
s
（11）
where λs (W/(m·°C)) is the thermal conductivity of formation. 2.2.5 Heat transfer model in formation Similarly, the energy balance equation for geothermal reservoir is established under cylindrical coordinates:
2T 1 T 2T s cs T 2 r r r 2 z s t
（12）
2.3 Model solution In Fig. 1, the solution domain is divided into subdomains, which are numbered in axial direction and radial directions. Then the finite difference method with full implicit scheme is adopted to discrete these differential equations. Subsequently, the Gauss-Seidel iterations is employed to solve the above equations, and MATLAB serves as the programming language. (1) The discrete equation for the model in central pipe is expressed as: n 1 n 1 n 1 n AT 1 1, j B1T1, j 1 C1T2, j D1T1, j
Where
C1
A1 B1 C1 D1 B1 ,
（13）
m qC pm z
C r2 1 D1 m pm 1 R , t
(2) The discrete equation for the model in annular space is expressed as:
A2T2,n j 1 B2T2,n j 11 C2T1,nj1 D2T3,n j 1 E2T2,n j Where A2 B2 C2 D2 E2 , B2
（14）
m qC pm z
mC pm r52 r42 1 C2 , D 2r5 h3 , E2 R 2 t (3) The discrete equation for the casing is expressed as:
A3T3,n j 1 B3T3,n j 11 C3T3,n j 11 D3T2,n j 1 E3T4,n j 1 F3T3,n j
（15）
ACCEPTED MANUSCRIPT
Where
A3 B3 C3 D3 E3 F3
C3
w r62 r52
z
F3
2
,
D3 2r5 h3
B3
6 r62 r52
z
,
E3 ,
2
2wc ln((r7 r6 ) / (r6 r5 ))
wCw r62 r52 t
(4) The discrete equation for the cement is expressed as:
A4T4,n j 1 B4T4,n j 11 C4T4,n j 11 D4T3,n j 1 E4T5,n j 1 F4T4,n j here A4 B4 C4 D4 E4 F4 , B4
C4
c r72 r62
z
2
, D4
（16）
c r72 r62
z
2
2wc ln((r7 r6 ) / (r6 r5 ))
c Cc r82 r72 2cs , F8 E8 t ln((r8 r7 ) / (r7 r6 ))
(5) The discrete equation for the formation is expressed as: n 1 n 1 n 1 n 1 n 1 n AT i i , j BiTi , j 1 CiTi , j 1 DiTi 1, j EiTi 1, j FT i i, j
Where Ai Bi Ci Di Ei Fi , Bi
Ci
Ei
i ri 2 ri 21
z
2
, Di
（17）
i ri 2 ri 21
z
2
2i 1,i ln((ri 3 ri 2 ) / (ri 2 ri 1 ))
2i ,i 1 ln((ri 4 ri 3 ) / (ri 3 ri 2 ))
, Fi
s Cs ri 2 ri 21 t
2.4 Boundary and initial conditions (1) A Neumann boundary is imposed at the ground surface, where the adiabatic condition is applied:
T z
0 z 0
（18）
(2) The boundary temperature of formation in radial direction is the original formation temperature:
TIn1，j T0 a j z
（19）
Where T0 is the ground surface temperature. a is the geothermal gradient. Δz is the axial grid length. j is the number of grid cells in axial direction. I is the quantity of all grids in radial direction. (3) The initial temperature of working fluid, central pipe and other medium is the original formation
ACCEPTED MANUSCRIPT temperature. 2.5 Model validation To determine the accuracy of mathematic model, the field data of HGP-A well in Hawai is adopted, and the temperature distribution of geothermal reservoir is shown in Fig. 2 [14, 23]. Note that there was a power failure during the test in Fig. 3, but its influence is only in a short time and the temperature variation trend remains basically unchanged. Consequently, the data obtained from the experiment is reasonable. It is clear the calculation results agree with measured data very well at each point, and the temperature difference at 7d is only 0.09 °C, indicating the model proposed in this paper is quite reliable.
Fig. 2. Temperature distribution of geothermal reservoir in Hawaii [14]
ACCEPTED MANUSCRIPT Fig. 3. Comparison between simulation result and experimental data [23] In this paper, the implicit finite-difference algorithm is adopted, and it is unconditionally stable. However, the calculation error is related to grid size and time step, and there is a need to determine the two parameters. Note that the absolute tolerance is set to 10-5, which means that when the error of all grid nodes is smaller than 10-5, the iteration would achieve convergence, and the result of model validation indicates it is reasonable. Fig. 4 shows that the outlet temperature has a greater change when the radial grid length is altered. The reason is that the thermal process from the formation to working fluid mainly occurs in radial direction, so the refined radial grid is needed to calculate the sharp change of temperature gradient around the borehole. Obviously, if radial grid length is less than 0.4 m, the result would close to its exact value, though the axial grid length ranges from 1m to 20 m. This means the mesh density has meet the accuracy requirement. Therefore, considering precision and computing time, the radial grid length and axial grid length are set to 0.2 m and 20 m, respectively. And the total number of grids eventually reach 62500. Fig. 5 depicts that the time step only affects the result within 128 hours, then the outlet temperatures under different time steps show a good consistency. However, if the time step exceeds 32 h, it would be difficult to analyze the value of the time points within a smaller time step. And this also brings a larger error, so time step of 4 h is adopted in this paper.
Fig. 4. Comparison of outlet temperature under different grid sizes.
ACCEPTED MANUSCRIPT
Fig. 5. Comparison of outlet temperature under different time steps 3. Model parameters In this paper, the research is based on geothermal field in Xiong'an New Area, Hebei Province, China, which aims at exploiting geothermal resources with a depth of 5000 m. And it is intended to employ CBHE for its cleanness and better efficiency. Due to high cost of deep drilling and well completion, there is a need to study the heat extraction performance of CBHE to optimize its design. Based on the Technical specification for geothermal well drilling released the People's Republic of China, CBHE is designed as the third spudding structure, and its specific parameters are shown in Table 1. And Table 2 lists the physical parameters of heat conduction medium. Based on the geological data of Xiong'an geothermal reservoir, the ground surface temperature is 20 °C, and the geothermal gradient is 0.03 °C/m. Additionally, due to the climatic conditions in the north of China, the heating time is usually from November 15th to March 15th of the next year. Hence, the computation time in our model is set to 4 months. Additionally, the thermophysical parameters of water vary with temperature, which has a noticeable influence on convective heat coefficient. Consequently, the convective heat coefficient of all grid cells must be updated at each time step. Table 3 presents relevant convection heat transfer coefficient under different temperature. And the CP and AS represents the central pipe and annular space, respectively. Table 1 Well structure parameters of CBHE Structure parameters
Surface casing
Intermediate casing
Production casing
Inner diameter(mm)
384.1
276.4
201.2
Outer diameter(mm)
406.4
298.4
219.1
depth(m)
480
1240
5000
Cement return height(m)
0
0
0
Table 2 Physical parameters of heat conduction medium Physical parameters
Casing
Conventional
Air (0.2 bar)
Formation
Central pipe
ACCEPTED MANUSCRIPT cement Density (kg/m3) Thermal capacity (J/(kg·°C)) Thermal conductivity (W/(m·°C))
8060
2140
0.225
2650
7850
400
2000
1013
920
465
43.75
0.7
0.026
2.8
43.5
Table 3 Convection heat transfer coefficient under different temperature Temperatur e T (°C)
Thermal conductivity λ ((W/(m·°C))
Prandtl number Pr
Reynolds
Nusselt
Convection heat
number
number
transfer coefficient
Re
Nu
h (W/(m2·°C))
CP
AS
CP
AS
CP
AS
25
0.6026
6.1721
195906.8
34169.74
1011.116
224.6641
12041.48
1205.543
50
0.6381
3.5897
314626.3
54876.61
1167.931
265.7456
14728.39
1509.994
75
0.6635
2.3941
449286.1
78363.74
1272.338
294.6343
16683.73
1740.782
100
0.6794
1.7459
595840.1
103925.5
1343.699
315.3067
18041.69
1907.563
125
0.6868
1.3649
746109.3
130135.2
1393.379
330.1486
18912.51
2019.111
150
0.6862
1.1383
885799.1
154499.6
1429.534
340.9457
19386.29
2083.321
4. Analysis of temperature field According to the optimization results in section 5, the flow rate and inlet temperature are set to 25 m3/h and 25 °C, respectively. Fig. 6 displays that the impact scope is limited by comparison with the size of formation in radial direction, which also indicates the computational domain adopted in our model is reasonable. Apparently, a much larger temperature gradient appears around the borehole, and that's why the thermal process from the reservoir to working fluid mainly occurs along radial direction. Besides, a local enlargement is carry out to distinguish the wellbore from the surrounding formation.
ACCEPTED MANUSCRIPT
Fig. 6. Temperature contours of the wellbore and surrounding formation at 120d According to different thermal insulation conditions, we investigate the temperature distribution of working fluid in wellbore at 120d. Fig. 7 shows that the temperature in annulus decreases slightly within a certain depth. This is because the inlet temperature is higher than the temperature in formation, leading to the heat transfer from working fluid to reservoir. Then it increases gradually with well depth, while the temperature in central tube only has a little reduction, which proves its well heat insulation effect. However, for the CBHE without thermal insulation, the temperature in wellbore changes sharply, which indicates that the heat transfer between the central pipe and annulus is much stronger than that near the wellbore. But this would cause much smaller temperature difference between the working fluid and formation, which is not beneficial for heat transfer. Thus, the difference of outlet temperature can even reach 32.71 °C. Additionally, if we employ partial insulation, the outlet temperature also has a remarkable increase in Fig. 7. Hence, further study is carried out in section 5 to analyze the feasibility of partial insulation for cutting costs of full thermal insulation.
ACCEPTED MANUSCRIPT
Fig. 7. Temperature distribution of working fluid in wellbore at 120d Because the impact scope is an especially important parameter to analyze well spacing and the heat extraction capacity of geothermal reservoir [33], then the temperature difference distribution with a depth of 4550 m in radial direction is evaluated. This can be obtained by subtraction between the initial temperature and the current temperature in formation. And further analysis is preformed to compare the temperature difference distribution with full thermal insulation and without thermal insulation. The impact scope represents the radial distance when the temperature difference is close to zero. Fig. 8 indicates that the impact scope with full heat insulation is larger than that without heat insulation, which can eventually reach 13 m after 120 days. And it is evident that the temperature gradient with full thermal insulation has a much more dramatic change near the wellbore, which is beneficial to heat extraction. Furthermore, as a contrast, the production time has little influence on the impact scope without thermal insulation, meaning its thermal process achieves stable faster. But this would lead to a lower heat output. Then the temperature difference distribution with full insulation along axial direction is analyzed in Fig. 9. It is clear the impact scope and temperature gradient increase with well depth. The reason is that a larger temperature difference between working fluid and the deep reservoir would enhance heat transfer.
ACCEPTED MANUSCRIPT
Fig. 8. Comparison of temperature difference distribution with full thermal insulation and without thermal insulation at 30d, 60d, 90d, 120d
Fig. 9. Temperature difference distribution with full thermal insulation at 1550 m, 2550 m, 3550 m, 4550 m To determine the impact scope and further analyze the sustainability of CBHE, the temperature difference distribution before the next space heating in buildings is analyzed in Fig. 10, just 8 months in length. Apparently, the temperature difference distribution in axial direction has the
ACCEPTED MANUSCRIPT similar trend changing with that in Fig. 9. Compared to the temperature difference at 4550 m after the first space heating, the temperature difference at 4550 m in Fig. 10 even decreases by 92.63 °C. This means that the temperature field in formation almost recovers its original condition, verifying the sustainability of CBHE. In addition, note that the impact scope has an increase of 14 m, which shows that the determination of the well spacing needs to consider the effect of recovery period. Hence, further study is carried out in Fig. 11. It is clear the impact scope increases with recovery period, which ranging from 21 m to 35 m. Furthermore, the decrease trend of temperature difference at 4550 m tends to be slow, and it only has maximum decline of 2.57 °C from 8th month to 20th month. This also indicates that there is no need to employ a longer recovery period for the next indoor heating.
Fig. 10. Temperature difference distribution at 1550 m, 2550 m, 3550 m, 4550 m after 8 months of recovery
ACCEPTED MANUSCRIPT
Fig. 11. Temperature difference distribution at 4550 m under different recovery periods 5. Analysis of heat extraction performance To achieve the best heat extraction performance, it is highly necessary to study the effects of the key factors, including inlet temperature, mass flow, insulation section length of central tube, thermal conductivity of cement and reservoir conditions, etc. And the formula for thermal power is as follows:
Qout 0.001q out Cout Tout in Cin Tin
（20）
5.1 The effect of flow rate Based on the ground surface condition, the inlet temperature is set to 25 °C. Fig. 12 illustrates that the heat output decreases drastically in a particularly short time because of the sharp change of temperature near the wellbore. Additionally, it is difficult to compensate the energy loss around the borehole in time. Hence, the decreasing rate of outlet temperature would slow down gradually, then remains relatively stable. And the outlet temperature declines with flow rates, which leads to a smaller difference between the outlet temperature and inlet temperature in Eq. 20. But due to the increase of the fluid mass participating in heat transfer in unit time, the thermal power would have an upward trend. In view of the heat supply requirement, the outlet temperature of CBHE needs to be kept above 60 °C to maintain the appropriate indoor temperature [34]. This means that the flow rate should be less than 25 m3/h in Fig. 13, and its outlet temperature is approximately 61.81 °C at 120d. Furthermore, the change flow rate has a dramatic effect on the pressure drop, which is a quite important parameter during production. Hence, the pressure drop under different flow rates is investigated in Fig. 14. And we employ the assumption of a smooth tube for the calculation, which has been proved its good precision for a U-tube collector [35]. It is clear the pressure drop has a noticeable increase with flow rate, which can range from 3.729 MPa to 22.021 MPa. Moreover, the growth rate of pressure drop increases as rises the flow rate, while that of thermal power slows down
ACCEPTED MANUSCRIPT gradually. Therefore, in terms of pressure drop and thermal power, 25 m3/h would be best for the CBHE performance. Due to the pressure drop at 25 m3/h satisfies the production requirements, so the flow rate of 25 m³/h is employed for the following analysis.
Fig. 12. Outlet temperature and thermal power with production time under different flow rates
Fig. 13. Outlet temperature and thermal power with flow rate at 30d, 60d, 90d, 120d
ACCEPTED MANUSCRIPT
Fig. 14. Thermal power and pressure drop with flow rate 5.2 The effect of inlet temperature Fig. 15. illustrates that the outlet temperature increases with inlet temperature, but the thermal power declines linearly, and its curve slope is approximately -13.55 kW/°C at 120d. This indicates that the inlet temperature has significant effect on the heat extraction performance of CBHE. And the lower temperature difference between the reservoir and working fluid would weaken the heat transfer. Accordingly, if the outlet temperature has meet the heating requirements, there is no need to raise the inlet temperature. Thus, the 25 °C serve as the inlet temperature in following study based on the surface condition.
Fig. 15. Outlet temperature and thermal power with inlet temperature at 30d, 60d, 90d, 120d 5.3 The effect of insulation section length Fig. 16 shows that the growth rate of outlet temperature and thermal power increases when the
ACCEPTED MANUSCRIPT insulation length is less than 3500 m, then this trend becomes slow gradually and remains basically unchanged after it exceeds 4500 m, verifying the feasibility of partial insulation. And its thermal power has only a drop of 1.12% compared to that with full insulation at 120d. The reason is that the temperature difference between annulus and central pipe declines gradually with well depth according to the temperature distribution of working fluid in Fig. 7, which can weaken the thermal process near the bottom. Additionally, the difference of outlet temperature and thermal power tends to be larger with insulation length. By comparison with the production curve in Fig. 12, it can be found that a shorter insulation length would make the production stable in less time and further lead to a drop of total energy output.
Fig. 16. Outlet temperature and thermal power with insulation section length at 30d, 60d, 90d,120d 5.4 The effect of cement sheath Because of the lower thermal conductivity of cement compared to that of formation, so it may hinder the heat transfer between the reservoir and working fluid. This phenomenon is very similar to the skin effect in petroleum engineering, which has significant influence on the oil production [36, 37]. Obviously, the outlet temperature and thermal power have a significant change when the thermal conductivity is less than 1.0 W/(m·K) in Fig. 17, then this trend slows down gradually, verifying the existence of skin effect in geothermal development. And the thermal power at 2.8 W/(m·K) (corresponding to the formation) can increase by 12.22% compared to that at 0.7 W/(m·K) after 120 days. Note that if the well cementing operations is not carried out, the formation water (0.569~0.687 W/(m·K) [38]) would fill its space, which can result in a worse performance. Apparently, the heat extraction region mainly lies in the deep formation, so there is no need to employ high thermal conductivity cement in whole length. Fig. 18 displays that the outlet temperature difference and thermal power difference with high thermal conductivity length. This is conducted by subtracting the outlet temperature and thermal power at a certain high thermal conductivity length from that with high thermal conductivity cement in whole length, respectively. The thermal conductivity of cement sheath is set to 2.8 W/(m·K). It should be emphasized that 0m on the x axis means the depth of well bottom, while 3760 m on the x axis represents the depth of
ACCEPTED MANUSCRIPT intermediate casing. It can be seen the trends of reduction of outlet temperature difference and thermal power difference tend to be slow with high thermal conductivity length. And there is only an outlet temperature difference of 0.31 °C at 3760 m after 120 days, which proves the feasibility of partial cementing with high thermal conductivity cement.
Fig. 17. Outlet temperature and thermal power with thermal conductivity of cement at 30d, 60d, 90d,120d
Fig. 18. Outlet temperature difference and thermal power difference with high thermal conductivity length at 30d, 60d, 90d,120d 5.5 The effect of geological conditions The geological conditions are particularly important basis for the location optimization of a geothermal well, so the effects of the thermal conductivity of formation and geothermal gradient is analyzed. Fig. 19 shows that the outlet temperature and thermal power increase with the thermal
ACCEPTED MANUSCRIPT conductivity of formation, and this trend slows down gradually. Furthermore, the thermal power at 4.0 W/(m·K) ((corresponding to granite[39]) is approximately 1290.24 kW compared to the 555.29 kW at 1.0 W/(m·K) (corresponding to shallow soil[39]) after 120 days, and its percentage of increase reaches 132.35%. This also indicates that the CBHE system is more suitable for the development of deep geothermal resources. Fig. 20 illustrates that the outlet temperature and thermal power increase linearly with geothermal gradient, and the curve slope can even reach 13.06 °C/(°C/100m) and 379.19 kW/(°C/100m) at 120d, respectively. Consequently, the thermal conductivity of formation and geothermal gradient both have remarkable impact on the heat extraction performance.
Fig. 19. Outlet temperature and thermal power with thermal conductivity of formation at 30d, 60d, 90d,120d
Fig. 20. Outlet temperature and thermal power with geothermal gradient at 30d, 60d, 90d,120d 5.6 The effect of well depth
ACCEPTED MANUSCRIPT For the cost assessment of a geothermal well, the well depth is a very important factor [40], so the thermal characteristics under different well depth are investigated to provide guidance for the well depth optimization. Fig. 21 depicts that the growth rate of outlet and thermal power increase with well depth because of the heat transfer enhancement between the reservoir and working fluid. And if we employ a well depth of 2000 m, the outlet temperature difference can reach 30.88 °C compared to the outlet temperature at 5000 m (corresponding the design depth). Hence, the CBHE system has great superiority in exploiting the deep geothermal resources. Furthermore, it would take longer to achieve stable as rises the well depth according to the production curve in Fig. 12, which can further lead to the increase of total energy output.
Fig. 21 Outlet temperature and thermal power with well depth at 30d, 60d, 90d,120d 6. Conclusion In this paper, an unsteady-state heat transfer model is proposed for a deep CBHE. The finite difference method is employed to solve the heat transfer model, and the model is validated with experiment data. Subsequently, the temperature field in formation is investigated. And the influences of key parameters on heat extraction performance of CBHE are analyzed. The main conclusions are as follows: (1) The temperature field analysis shows that there is a significant temperature gradient around the borehole, far more than the geothermal gradient, indicating the heat transfer mainly occurs in the radial direction. The impact scope in the reservoir is relatively limited, so the CBHE has a great potential in developing geothermal resources. Furthermore, the temperature field in formation before the next heating period would get close to its initial condition, which illustrates that the heat production of CBHE system is sustainable. And the impact scope increases with recovery period, which is an important reference for well spacing optimization. (2) The outlet temperature and thermal power have a remarkable decrease at the initial stage, but then remains relatively stable. Moreover, a higher flow rate leads to the decline of outlet temperature, but the thermal power would increase. Further study indicates that in terms of the pressure drop and thermal power, the 25 m3/s is the optimal flow rate which could obtain the best CBHE performance. (3) As the inlet temperature increases, the outlet temperature declines and the thermal power increases. This depicts that a larger temperature difference between the reservoir and working fluid
ACCEPTED MANUSCRIPT could enhance heat transfer. Therefore, a lower inlet temperature is recommended if the outlet temperature has meet the required temperature. (4) The temperature of working fluid in central tube only has a little reduction, verifying its well heat insulation effect. Additionally, the insulation section length of central tube only has a significant effect on heat extraction performance within a certain range. And the changes of outlet temperature and thermal power are minor when this length exceeds 4500 m, which indicates that there is an optimal insulation length. (5) A lower thermal conductivity of cement sheath could weaken the heat transfer from the reservoir to working fluid, verifying the existence of the skin effect in geothermal development. Consequently, it is necessary to conduct the research of high thermal conductivity cement. And to reduces its cost, further analysis depicts that there is an optimum length with high thermal conductivity cement because the thermal process mainly lies in the deep reservoir. (6) The increase of the thermal conductivity of formation and geothermal gradient can greatly enhance the heat exchange, so it is particularly important to determine the appropriate geological conditions for CBHE system. Additionally, the well depth has a noticeable influence on the thermal process, and even the growth rate of outlet temperature and thermal power increase with well depth. The findings can offer guidance for the location optimization and well depth optimization of CBHE system. Acknowledgements The authors would like to acknowledge the National Key Research and Development Program of China (Grant No. 2016YFE0124600). Besides, support from the Program of Introducing Talents of Discipline to Chinese Universities (111 Plan) (Grant NO. B17045) is appreciated. Nomenclature AS = annular space a = geothermal gradient, °C/m CBHE = coaxial borehole heat exchanger CP = central pipe Cpm = heat capacity of working fluid, J/(kg·°C) Cin = heat capacity of inlet fluid, J/(kg·°C) Cout = heat capacity of outlet fluid, J/(kg·°C) Cw = heat capacity of casing, J/(kg·°C) Cc = heat capacity of cement sheath, J/(kg·°C) Cs = heat capacity of formation, J/(kg·°C) Deq = hydraulic diameter, m f = Darcy friction factor h1 = convective heat transfer coefficient for inner surface of central pipe, W/(m2·°C) h2 = convective heat transfer coefficient for outer surface of central pipe, W/(m2·°C) h3 = convective heat transfer coefficient for casing wall, W/(m2·°C) Nu = Nusselt number Pr = Prandtl number Qout = thermal power, kW q = flow rate of working fluid, m³/h R = thermal resistance, (m·°C)/W Re = Reynolds number
ACCEPTED MANUSCRIPT r1 r2 r3 r4 r5 r6 r7 Tf1 Tf2 Tin Tout Tw Tc ρm ρin ρout ρw ρc ρs λf λw1 λw2 λw3 λw λc λwc λs λcs Δt Δz i j
= = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
inner radius of central pipe, m inner radius of insulation layer, m outer radius of insulation layer, m outer radius of central pipe, m inner radius of casing, m outer radius of casing, m outer radius of cement sheath, m temperature of working fluid in central pipe, °C temperature of working fluid in annulus space, °C temperature of inlet fluid, °C temperature of outlet fluid, °C temperature of casing, °C temperature of cement sheath, °C density of working fluid, kg/m3 density of inlet fluid, kg/m3 density of outlet fluid, kg/m3 density of casing, kg/m3 density of cement sheath, kg/m3 density of formation, kg/m3 thermal conductivity of working fluid, W/(m·°C) thermal conductivity of inner layer of central pipe, W/(m·°C) thermal conductivity of insulation layer of central pipe, W/(m·°C) thermal conductivity of outer layer of central pipe, W/(m·°C) thermal conductivity of casing, W/(m·°C) thermal conductivity of cement sheath, W/(m·°C) equivalent thermal conductivity of casing and cement, W/(m·°C) thermal conductivity of formation, W/(m·°C) equivalent thermal conductivity of cement and formation, W/(m·°C) time step, s length of axial grid, m number of grid cells in radial direction number of grid cells in axial direction
References [1] J.C. Reboredo, Renewable energy contribution to the energy supply: Is there convergence across countries?, Renewable & Sustainable Energy Reviews 45 (2015) 290-295. [2] M. Sun, P.P. Zhang, T.H. Shan, C.C. Fang, X.F. Wang, L.X. Tian, Research on the evolution model of an energy supply– demand network, Physica A Statistical Mechanics & Its Applications 391(19) (2012) 4506-4516. [3] M.E. Bildirici, S.M. Gökmenoğlu, Environmental pollution, hydropower energy consumption and economic growth: Evidence from G7 countries, Renewable & Sustainable Energy Reviews 75 (2017) 68-85. [4] Z. Lyu, X. Song, G. Li, X. Hu, Y. Shi, Z. Xu, Numerical Analysis of Characteristics of a Single U-tube Downhole Heat Exchanger in the Borehole for Geothermal Wells, Energy 125 (2017) 186-196. [5] R.A. Beier, Transient heat transfer in a U-tube borehole heat exchanger, Applied Thermal Engineering 62(1) (2014) 256-266. [6] H. Holmberg, J. Acuña, E. Næss, O.K. Sønju, Thermal evaluation of coaxial deep borehole heat exchangers, Renewable Energy
ACCEPTED MANUSCRIPT 97 (2016) 65-76. [7] A. Zarrella, M. Scarpa, M.D. Carli, Short time-step performances of coaxial and double U-tube borehole heat exchangers: Modeling and measurements, Hvac & R Research 17(6) (2011) 959-976. [8] R.N. Horne, Design considerations of a down-hole coaxial geothermal heat exchanger, Trans. - Geotherm. Resour. Counc.; (United States) 4 (1980). [9] P.J. Yekoladio, T. Bello-Ochende, J.P. Meyer, Design and optimization of a downhole coaxial heat exchanger for an enhanced geothermal system (EGS), Renewable Energy 55(4) (2013) 128-137. [10] R.A. Beier, J. Acuña, P. Mogensen, B. Palm, Transient heat transfer in a coaxial borehole heat exchanger, Geothermics 51(7) (2014) 470-482. [11] D. Gordon, T. Bolisetti, S.K. Ting, S. Reitsma, Short-term fluid temperature variations in either a coaxial or U-tube borehole heat exchanger, Geothermics 67 (2017) 29-39. [12] D. Gordon, T. Bolisetti, S.K. Ting, S. Reitsma, A physical and semi-analytical comparison between coaxial BHE designs considering various piping materials, Energy 141 (2018). [13] K. Morita, T. Yamaguchi, H. Karasawa, H. Hayamizu, Development and evaluation of temperature simulation code for geothermal wells, J. Min. Metall. Inst. Jpn 100 (1984) 1045-1051. [14] K. Morita, W.S. Bollmeier, H. Mizogami, Analysis of the results from the Downhole Coaxial Heat Exchanger (DCHE) experiment in Hawaii,
(1992).
[15] M.S. Saadi, R. Gomri, Investigation of dynamic heat transfer process through coaxial heat exchangers in the ground, International Journal of Hydrogen Energy 42(28) (2017). [16] E. Zanchini, S. Lazzari, A. Priarone, Improving the thermal performance of coaxial borehole heat exchangers, Energy 35(2) (2010) 657-666. [17] E. Zanchini, S. Lazzari, A. Priarone, Effects of flow direction and thermal short-circuiting on the performance of small coaxial ground heat exchangers, Renewable Energy 35(6) (2010) 1255-1265. [18] X. Song, Y. Shi, G. Li, R. Yang, Z. Xu, R. Zheng, G. Wang, Z. Lyu, Heat Extraction Performance Simulation for Various Configurations of a Downhole Heat Exchanger Geothermal System, Energy 141 (2017). [19] Y. Shi, X. Song, G. Li, R. Li, Y. Zhang, G. Wang, R. Zheng, Z. Lyu, Numerical Investigation on Heat Extraction Performance of a Downhole Heat Exchanger Geothermal System, Applied Thermal Engineering
(2018).
[20] Y. Fujimitsu, K. Fukuoka, S. Ehara, H. Takeshita, F. Abe, Evaluation of subsurface thermal environmental change caused by a ground-coupled heat pump system, Current Applied Physics 10(2) (2010) S113-S116. [21] K. Bär, W. Rühaak, B. Welsch, D. Schulte, S. Homuth, I. Sass, Seasonal High Temperature Heat Storage with Medium Deep Borehole Heat Exchangers ☆, Energy Procedia 76 (2015) 351-360. [22] G. Cui, S. Ren, L. Zhang, J. Ezekiel, C. Enechukwu, Y. Wang, R. Zhang, Geothermal exploitation from hot dry rocks via recycling heat transmission fluid in a horizontal well, Energy 128 (2017) 366-377. [23] Morita, I I u AN EXPERIMENT TO PROVE THE CONCEPT OF THE DOWNHOLE COAXIAL HEAT EXCHANGER (DCHE) IN HAWAII, Geothermal Resources Council Trans 16 (1992) 9-16. [24] T. Kohl, M. Salton, L. Rybach, Data analysis of the deep borehole heat exchanger plant weissbad (Switzerland), Proc World Geothermal Congress
(2000).
[25] T. Kohl, R. Brenni, W. Eugster, System performance of a deep borehole heat exchanger, Geothermics 31(6) (2002) 687-708. [26] J. Acuña, B. Palm, Distributed thermal response tests on pipe-in-pipe borehole heat exchangers, Applied Energy 109(6) (2013) 312-320. [27] J. Zhu, K. Hu, X. Lu, X. Huang, K. Liu, X. Wu, A review of geothermal energy resources, development, and applications in China: Current status and prospects, Energy 93 (2015) 466-483. [28] J.W. Lund, T.L. Boyd, Direct utilization of geothermal energy 2015 worldwide review, Geothermics 60 (2016) 66-93. [29] J. Hou, M. Cao, P. Liu, Development and utilization of geothermal energy in China: Current practices and future strategies,
ACCEPTED MANUSCRIPT Renewable Energy 125 (2018). [30] A.D. Chiasson, G.G. Culver, D. Favata, S. Keiffer, DESIGN, INSTALLATION, AND MONITORING OF A NEW DOWNHOLE HEAT EXCHANGER, IEEE International Conference on Granular Computing, 2005. [31] V. Gnielinski, New equations for heat and mass transfer in the turbulent flow in pipes and channels, Nasa Sti/recon Technical Report A 75(2) (1975) 8-16. [32] G.K. Filonenko, Hydraulic Resistance in Pipes, Teploenergetika 4(4) (1954) 15-21. [33] Y. Shi, X. Song, G. Li, Z. Shen, X. Hu, Z. Lyu, R. Zheng, G. Wang, Numerical Analysis of the Heat Production Performance of a Closed Loop Geothermal System, Renewable Energy 120 (2017). [34] T. Nussbaumer, S. Thalmann, Status report on district heating systems in IEA countries, Swiss Federal Office of Energy, and Verenum, Zürich
(2014).
[35] J. Acuña, Improvements of U-pipe Borehole Heat Exchangers, Annan Naturresursteknik (2010). [36] A.F.V. Everdingen, The Skin Effect and Its Influence on the Productive Capacity of a Well, Journal of Petroleum Technology 5(6) (1953) 171-176. [37] W. Hurst, J.D. Clark, E.B. Brauer, The skin effect in producing wells, Journal of the Society of Petroleum Evaluation Engineers (1968). [38] E.O. Holzbecher, Modeling Density-Driven Flow in Porous Media, Springer1998. [39] G. Hellström, Ground heat storage: Thermal analyses of duct storage systems, Physics & Astronomy
(1991).
[40] LUKAWSKI M. Z, ANDERSON B. J, AUGUSTINE Chad, CAPUANO L. E, BECKERS K. F, LIVESAY Bill, TESTER J. W, Cost Analysis of Oil, Gas, and Geothermal Well Drilling, Journal of Petroleum Science & Engineering 118 (2014) 1-14.
ACCEPTED MANUSCRIPT
Highlights of “Numerical analysis of heat extraction performance of a deep coaxial borehole heat exchanger geothermal system” Xianzhi Song1*, Gaosheng Wang1, Yu Shi1, Ruixia Li2, Zhengming Xu1, Rui Zheng1, Yu Wang1, Jiacheng Li1 (1. State Key Laboratory of Petroleum Resources and Prospecting, China University of Petroleum, Beijing, Beijing 102249, China; 2. Sinopec Star Petroleum Co. Ltd, Beijing 100083, China)
An unsteady-state heat transfer numerical model for CBHE is presented. The simulation values show particularly good agreement with field data. The temperature field dynamic response of CBHE is analyzed comprehensively. Influences of key factors on the heat production performance of CBHE are studied.
Xianzhi Song, Gaosheng Wang, Yu Shi, Gensheng Li, Ruixia Li, Zhengming Xu, Rui Zheng, Yu Wang PII:
S0360-5442(18)31588-3
DOI:
10.1016/j.energy.2018.08.056
Reference:
EGY 13532
To appear in:
Energy
Received Date:
09 April 2018
Accepted Date:
07 August 2018
Please cite this article as: Xianzhi Song, Gaosheng Wang, Yu Shi, Gensheng Li, Ruixia Li, Zhengming Xu, Rui Zheng, Yu Wang, Numerical Analysis of Heat Extraction Performance of a Deep Coaxial Borehole Heat Exchanger Geothermal System, Energy (2018), doi: 10.1016/j.energy. 2018.08.056
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT Numerical Analysis of Heat Extraction Performance of a Deep Coaxial Borehole Heat Exchanger Geothermal System Xianzhi Song1*, Gaosheng Wang1, Yu Shi1, Gensheng Li1, Ruixia Li2, Zhengming Xu1, Rui Zheng1, Yu Wang1 (1. State Key Laboratory of Petroleum Resources and Prospecting, China University of Petroleum, Beijing, Beijing 102249, China; 2. Sinopec Star Petroleum Co. Ltd, Beijing 100083, China)
ABSTRACT: The coaxial borehole heat exchanger (CBHE) can effectively develop deep geothermal resources with larger heat exchange area than traditional U-tube BHE. However, few studies focus on its parameters optimization and temperature field dynamic response. Hence, an unsteady-state heat transfer model is established for a deep CBHE. The finite difference method is adopted to solve the model, and the model is validated with experiment data. The effects of key factors on the CBHE performance are comprehensively analyzed. Simulation results depict that the outlet temperature decreases greatly at the initial stage, and then remains relatively stable. There is a critical value for flow rate to obtain higher thermal power with appropriate pressure drop. The insulated pipe can significantly reduce heat dissipation of working fluid in central pipe. There is an optimal insulation length in terms of heat insulation effect and costs. The thermal conductivity of cement sheath has remarkable influence on the thermal process and partial cementing with high thermal conductivity cement is favorable. The temperature field analysis shows that the impact scope is limited, and the temperature in the reservoir before the next heating period would almost recover its original condition, verifying the sustainability of CBHE system. KEY WORDS: geothermal energy; coaxial borehole heat exchanger; unsteady-state heat transfer; numerical simulation, heat extraction performance 1. Introduction The massive consume of the fossil fuels have caused many problems, such as resource depletion, global warming and environmental pollution [1, 2]. Therefore, the renewable energy such as geothermal energy is playing an increasingly important role in the district energy supply and environmental improvement [3]. The borehole heat exchanger has been widely used to extract geothermal energy for space heating in the residential and commercial buildings, and the U-tube borehole heat exchanger is the most common form [4, 5]. However, the coaxial borehole heat exchanger (CBHE) is preferred to develop the deep geothermal resources for larger heat exchange area and better performance [6, 7]. Investigation on the heat transfer process of CBHE system has been carried out by many researchers. Roland N. Horne [8] first established a quasi-steady heat transfer analytical model, and studied the effects of flow direction, inlet temperature and outer tube radius. To optimize the flow rate and geometric size, P.J. Yekoladio et al. [9] developed a steady state analytical model, while the dynamic change of temperature field in formation was ignored. Richard A. Beier et al. [10] proposed a transient heat transfer model solved by Laplace transform to investigate the borehole thermal resistance. David Gordon et al. [11, 12] presented composite cylindrical source (CCS) model to analyze the short-term performance of CBHE system, then the influence of pipe materials was evaluated by a semi-analytical model. Analytical models contribute to explore the rigorous relationship of each variable, while the numerical solution is also a convenient tool to study the performance of CBHE. Morita et al. [13, 14] employed explicit finite-difference method to solve the heat transfer model, and investigated the thermal characteristics. Henrik Holmberg et al. [6] discretized the energy equation by implicit finite difference scheme, and results showed that the * Corresponding author. E-mail address: [email protected] (X. Song).
ACCEPTED MANUSCRIPT deeper boreholes was beneficial for heat extraction. Then Mohamed Salah Saadi et al. [15] built a transient 2-D FVM (Finite volume method) model considering the seasonal effects for CBHE. Previous studies have made significant contributions to comprehend the thermal process of CBHE. However, few studies focus on the parameters optimization, for example the critical values of flow rate, thermal insulation length, cementing length, etc. And the investigation on the temperature field dynamic response in reservoir is still lacking, especially after a recovery period, which is vital for well spacing determination and sustainability analysis. Additionally, the commercial simulation software, such as COMSOL and FEFLOW, also provides a new platform to study the thermal process. Zanchini et al. [16, 17] established a 2-D model through COMSOL, then found that the thermal conductivity of cement played a remarkable role in the heat transfer. Xianzhi Song et al. [18, 19] presented a 3-D transient state model coupled with groundwater flow by COMSOL, and evaluated the heat transfer characteristics of different borehole heat exchangers. Based on FEFLOW code, some researchers analyzed the subsurface thermal environmental changes and seasonal heat storage efficiency [20, 21]. Nevertheless, the mesh generation for many thin layers of large scale model, such as casing and cement sheath, is not convenient in commercial software. Besides, Guodong Cui et al. [22] proposed a promising CBHE with single horizontal well, but the vertical CBHE is recommended in view of drilling cost and complexity. To evaluate the actual heat extraction performance of CBHE, a field experiment was carried out at HGP-A well in Hawaii, which also provided reliable data for model validation [14, 23]. In 1993, a deep CBHE in Weissbad, Switzerland, was operated commercially, but its outlet temperature was lower than expected. T. Kohl1 et al. [24, 25] attributed this to the inner tube without insulation and the poor contact between casing and formation. Furthermore, a distributed thermal response test (DTRT) for CBHE was conducted to accurately analyze borehole thermal resistances [26]. China has the most installed thermal capacity for geothermal direct-use, but its geothermal power development is far below than that of the top 10 countries, such as United States, Mexico and Iceland [27-29]. The main reason is the relatively low temperature in single geothermal well, which leads to a limited heat output. Thus, there is an urgent need to exploit the deep geothermal resources in China. And the CBHE with a vertical wellbore would be employed to extract geothermal energy in Xiong'an New Area for its better performance and cleanliness. However, there are few comprehensive researches of parameters optimization, including flow rate, inlet temperature, cementing length, insulation section length, especially for a deep CBHE. And it is still lack of further analysis on the dynamic change of temperature field in formation. Hence, an unsteady-state heat transfer model is proposed in this paper. Then the finite difference method is adopted to solve the mathematical model, and the model is validated with experiment data. The main findings can offer guidance for the optimal design of CBHE system. 2. Model establishment 2.1 Physical Model The CBHE system consists of high-pressure pump, vertical wellbore and geothermal heat pump in Fig. 1. A working fluid is injected from the annulus, and then produced through central pipe. Forced convection occurs in the wellbore, while the thermal process in solid simply involves heat conduction. Through heat pump, the geothermal energy carried by the high-temperature fluid can be used for space heating or power generation. Additionally, because of the high-pressure environment in a deep well, the conventional BHE materials with low thermal conductivity [6, 30],
ACCEPTED MANUSCRIPT such as polyethylene and polypropylene, are difficult to meet the strength requirement. Hence, a three-layer vacuum tube is employed, and the pressure of insulation layer remains 0.2 bar by a vacuum pump, which avoids air circulation.
Fig. 1. The schematic diagram of the flow and heat transfer in CBHE 2.2 Mathematic model Due to the symmetry of simulation domain and boundary conditions, a 2D unsteady-state heat transfer model is set up based on the energy balance relationship among the heat exchange medium. And the circulating fluid, inner tube, casing and cement sheath are simplified into one-dimension, while we adopt two-dimensional grid for the reservoir simulation. And the assumptions in this model include: (1) The heat generation by viscous friction of working fluid is negligible compared to the heat exchange between working fluid and its surrounding medium; (2) Regarding water as the working fluid, and its heat transfer form is limited to forced convection; (3) Only conductive heat transport is considered in formation. 2.2.1 Heat transfer model in central tube According to Fig. 1, the energy balance equation is established for the working fluid. The second term on the left side refers to heat transfer from annular fluid. The thermal process involves forced convection and heat conduction through central tube. The expression is as follows:
m qC pm
T f 1 z
Tf 2 Tf 1 R
mC pm r12
T f 1 t
（1）
where Tf1 and Tf2 (°C) are the temperature of working fluid in central pipe and annulus, respectively. R ((m·°C)/W) is the thermal resistance, which is similar as electric resistance. This variable is defined as:
R (W/(m2·°C))
ln r2 / r1 ln r3 / r2 ln r4 / r3 1 1 2 r1h1 2w1 2w 2 2w3 2 r4 h2
（2）
where h is the convective heat transfer coefficient, which is affected by temperature due to the change of physical properties of working fluid. This coefficient is defined as:
ACCEPTED MANUSCRIPT
h
Nu f D
（3）
where Nu is Nusselt number, which characterizes the intensity of convective heat transfer. It is calculated by Gnielinski equation [31]:
f Re 1000 Pr 8 Nu 1 12.7 f
8
Pr
23
1
（4） where Re is Reynolds number. Pr is Prandtl number. f is Darcy friction factor of turbulent flow in a tube, which is calculated by Filonenko equation [32]:
f
1
1.82 lg Re 1.64
2
（5）
Re=3000~6×106,
The application range of Eq. (5) is: Pr=0.5~2000, L/D>10, which represents the ratio of casing length to diameter. Then the hydraulic diameter is introduced for the flow in annulus:
Deq 2r5 2r4
（6）
where r4 and r5 (m) are the inner and outer radius of the annular space, respectively. 2.2.2 Heat transfer model in annular space The energy balance equation for annular fluid is established based on the thermal process. Compared with Eq. (1), there is one more item on the left side in Eq. (7), which represents the heat transfer from casing, as shown below:
m qC pm
T f 2 z
Tf 1 Tf 2 R
2 r5 h3 Tw T f 2 mC pm r52 r42
T f 2 t
（7）
2.2.3 Heat transfer model in casing Due to the temperature difference, the heat would transfer from cement sheath to casing, which is represented by the third item in Eq. (8). Also, it needs to consider the forced convection heat transfer. The equation is described by:
2Tw 2wc w r r 2 2 r5 h3 T f 2 Tw Tc Tw z ln((r7 r6 ) / (r6 r5 )) 2 6
2 5
wCw r62 r52
Tw t
（8）
where λwc (W/(m·°C)) is the equivalent thermal conductivity of casing and cement, which is determined by:
wc =
ln (r7 r6 ) (r6 r5 ) ln 2r6 (r6 r5 ) ln (r7 r6 ) 2r6
w
c
（9）
where λw and λc (W/(m·°C)) are the thermal conductivity of casing and cement sheath, respectively. 2.2.4 Heat transfer model in cement sheath The heat transfer process in cement sheath is only limited to heat conduction, and its energy balance
ACCEPTED MANUSCRIPT equation is set up as:
2Tc 2wc 2cs c r r 2 Tw Tc Ts Tc z ln((r7 r6 ) / (r6 r5 )) ln((r8 r7 ) / (r7 r6 )) 2 7
2 6
c Cc r72 r62
Tc t
（10）
where λcs (W/(m·°C)) is the equivalent thermal conductivity of cement and formation, and it is identified by:
cs =
ln (r8 r7 ) (r7 r6 ) ln 2r7 (r7 r6 ) ln (r8 r7 ) 2r7
c
s
（11）
where λs (W/(m·°C)) is the thermal conductivity of formation. 2.2.5 Heat transfer model in formation Similarly, the energy balance equation for geothermal reservoir is established under cylindrical coordinates:
2T 1 T 2T s cs T 2 r r r 2 z s t
（12）
2.3 Model solution In Fig. 1, the solution domain is divided into subdomains, which are numbered in axial direction and radial directions. Then the finite difference method with full implicit scheme is adopted to discrete these differential equations. Subsequently, the Gauss-Seidel iterations is employed to solve the above equations, and MATLAB serves as the programming language. (1) The discrete equation for the model in central pipe is expressed as: n 1 n 1 n 1 n AT 1 1, j B1T1, j 1 C1T2, j D1T1, j
Where
C1
A1 B1 C1 D1 B1 ,
（13）
m qC pm z
C r2 1 D1 m pm 1 R , t
(2) The discrete equation for the model in annular space is expressed as:
A2T2,n j 1 B2T2,n j 11 C2T1,nj1 D2T3,n j 1 E2T2,n j Where A2 B2 C2 D2 E2 , B2
（14）
m qC pm z
mC pm r52 r42 1 C2 , D 2r5 h3 , E2 R 2 t (3) The discrete equation for the casing is expressed as:
A3T3,n j 1 B3T3,n j 11 C3T3,n j 11 D3T2,n j 1 E3T4,n j 1 F3T3,n j
（15）
ACCEPTED MANUSCRIPT
Where
A3 B3 C3 D3 E3 F3
C3
w r62 r52
z
F3
2
,
D3 2r5 h3
B3
6 r62 r52
z
,
E3 ,
2
2wc ln((r7 r6 ) / (r6 r5 ))
wCw r62 r52 t
(4) The discrete equation for the cement is expressed as:
A4T4,n j 1 B4T4,n j 11 C4T4,n j 11 D4T3,n j 1 E4T5,n j 1 F4T4,n j here A4 B4 C4 D4 E4 F4 , B4
C4
c r72 r62
z
2
, D4
（16）
c r72 r62
z
2
2wc ln((r7 r6 ) / (r6 r5 ))
c Cc r82 r72 2cs , F8 E8 t ln((r8 r7 ) / (r7 r6 ))
(5) The discrete equation for the formation is expressed as: n 1 n 1 n 1 n 1 n 1 n AT i i , j BiTi , j 1 CiTi , j 1 DiTi 1, j EiTi 1, j FT i i, j
Where Ai Bi Ci Di Ei Fi , Bi
Ci
Ei
i ri 2 ri 21
z
2
, Di
（17）
i ri 2 ri 21
z
2
2i 1,i ln((ri 3 ri 2 ) / (ri 2 ri 1 ))
2i ,i 1 ln((ri 4 ri 3 ) / (ri 3 ri 2 ))
, Fi
s Cs ri 2 ri 21 t
2.4 Boundary and initial conditions (1) A Neumann boundary is imposed at the ground surface, where the adiabatic condition is applied:
T z
0 z 0
（18）
(2) The boundary temperature of formation in radial direction is the original formation temperature:
TIn1，j T0 a j z
（19）
Where T0 is the ground surface temperature. a is the geothermal gradient. Δz is the axial grid length. j is the number of grid cells in axial direction. I is the quantity of all grids in radial direction. (3) The initial temperature of working fluid, central pipe and other medium is the original formation
ACCEPTED MANUSCRIPT temperature. 2.5 Model validation To determine the accuracy of mathematic model, the field data of HGP-A well in Hawai is adopted, and the temperature distribution of geothermal reservoir is shown in Fig. 2 [14, 23]. Note that there was a power failure during the test in Fig. 3, but its influence is only in a short time and the temperature variation trend remains basically unchanged. Consequently, the data obtained from the experiment is reasonable. It is clear the calculation results agree with measured data very well at each point, and the temperature difference at 7d is only 0.09 °C, indicating the model proposed in this paper is quite reliable.
Fig. 2. Temperature distribution of geothermal reservoir in Hawaii [14]
ACCEPTED MANUSCRIPT Fig. 3. Comparison between simulation result and experimental data [23] In this paper, the implicit finite-difference algorithm is adopted, and it is unconditionally stable. However, the calculation error is related to grid size and time step, and there is a need to determine the two parameters. Note that the absolute tolerance is set to 10-5, which means that when the error of all grid nodes is smaller than 10-5, the iteration would achieve convergence, and the result of model validation indicates it is reasonable. Fig. 4 shows that the outlet temperature has a greater change when the radial grid length is altered. The reason is that the thermal process from the formation to working fluid mainly occurs in radial direction, so the refined radial grid is needed to calculate the sharp change of temperature gradient around the borehole. Obviously, if radial grid length is less than 0.4 m, the result would close to its exact value, though the axial grid length ranges from 1m to 20 m. This means the mesh density has meet the accuracy requirement. Therefore, considering precision and computing time, the radial grid length and axial grid length are set to 0.2 m and 20 m, respectively. And the total number of grids eventually reach 62500. Fig. 5 depicts that the time step only affects the result within 128 hours, then the outlet temperatures under different time steps show a good consistency. However, if the time step exceeds 32 h, it would be difficult to analyze the value of the time points within a smaller time step. And this also brings a larger error, so time step of 4 h is adopted in this paper.
Fig. 4. Comparison of outlet temperature under different grid sizes.
ACCEPTED MANUSCRIPT
Fig. 5. Comparison of outlet temperature under different time steps 3. Model parameters In this paper, the research is based on geothermal field in Xiong'an New Area, Hebei Province, China, which aims at exploiting geothermal resources with a depth of 5000 m. And it is intended to employ CBHE for its cleanness and better efficiency. Due to high cost of deep drilling and well completion, there is a need to study the heat extraction performance of CBHE to optimize its design. Based on the Technical specification for geothermal well drilling released the People's Republic of China, CBHE is designed as the third spudding structure, and its specific parameters are shown in Table 1. And Table 2 lists the physical parameters of heat conduction medium. Based on the geological data of Xiong'an geothermal reservoir, the ground surface temperature is 20 °C, and the geothermal gradient is 0.03 °C/m. Additionally, due to the climatic conditions in the north of China, the heating time is usually from November 15th to March 15th of the next year. Hence, the computation time in our model is set to 4 months. Additionally, the thermophysical parameters of water vary with temperature, which has a noticeable influence on convective heat coefficient. Consequently, the convective heat coefficient of all grid cells must be updated at each time step. Table 3 presents relevant convection heat transfer coefficient under different temperature. And the CP and AS represents the central pipe and annular space, respectively. Table 1 Well structure parameters of CBHE Structure parameters
Surface casing
Intermediate casing
Production casing
Inner diameter(mm)
384.1
276.4
201.2
Outer diameter(mm)
406.4
298.4
219.1
depth(m)
480
1240
5000
Cement return height(m)
0
0
0
Table 2 Physical parameters of heat conduction medium Physical parameters
Casing
Conventional
Air (0.2 bar)
Formation
Central pipe
ACCEPTED MANUSCRIPT cement Density (kg/m3) Thermal capacity (J/(kg·°C)) Thermal conductivity (W/(m·°C))
8060
2140
0.225
2650
7850
400
2000
1013
920
465
43.75
0.7
0.026
2.8
43.5
Table 3 Convection heat transfer coefficient under different temperature Temperatur e T (°C)
Thermal conductivity λ ((W/(m·°C))
Prandtl number Pr
Reynolds
Nusselt
Convection heat
number
number
transfer coefficient
Re
Nu
h (W/(m2·°C))
CP
AS
CP
AS
CP
AS
25
0.6026
6.1721
195906.8
34169.74
1011.116
224.6641
12041.48
1205.543
50
0.6381
3.5897
314626.3
54876.61
1167.931
265.7456
14728.39
1509.994
75
0.6635
2.3941
449286.1
78363.74
1272.338
294.6343
16683.73
1740.782
100
0.6794
1.7459
595840.1
103925.5
1343.699
315.3067
18041.69
1907.563
125
0.6868
1.3649
746109.3
130135.2
1393.379
330.1486
18912.51
2019.111
150
0.6862
1.1383
885799.1
154499.6
1429.534
340.9457
19386.29
2083.321
4. Analysis of temperature field According to the optimization results in section 5, the flow rate and inlet temperature are set to 25 m3/h and 25 °C, respectively. Fig. 6 displays that the impact scope is limited by comparison with the size of formation in radial direction, which also indicates the computational domain adopted in our model is reasonable. Apparently, a much larger temperature gradient appears around the borehole, and that's why the thermal process from the reservoir to working fluid mainly occurs along radial direction. Besides, a local enlargement is carry out to distinguish the wellbore from the surrounding formation.
ACCEPTED MANUSCRIPT
Fig. 6. Temperature contours of the wellbore and surrounding formation at 120d According to different thermal insulation conditions, we investigate the temperature distribution of working fluid in wellbore at 120d. Fig. 7 shows that the temperature in annulus decreases slightly within a certain depth. This is because the inlet temperature is higher than the temperature in formation, leading to the heat transfer from working fluid to reservoir. Then it increases gradually with well depth, while the temperature in central tube only has a little reduction, which proves its well heat insulation effect. However, for the CBHE without thermal insulation, the temperature in wellbore changes sharply, which indicates that the heat transfer between the central pipe and annulus is much stronger than that near the wellbore. But this would cause much smaller temperature difference between the working fluid and formation, which is not beneficial for heat transfer. Thus, the difference of outlet temperature can even reach 32.71 °C. Additionally, if we employ partial insulation, the outlet temperature also has a remarkable increase in Fig. 7. Hence, further study is carried out in section 5 to analyze the feasibility of partial insulation for cutting costs of full thermal insulation.
ACCEPTED MANUSCRIPT
Fig. 7. Temperature distribution of working fluid in wellbore at 120d Because the impact scope is an especially important parameter to analyze well spacing and the heat extraction capacity of geothermal reservoir [33], then the temperature difference distribution with a depth of 4550 m in radial direction is evaluated. This can be obtained by subtraction between the initial temperature and the current temperature in formation. And further analysis is preformed to compare the temperature difference distribution with full thermal insulation and without thermal insulation. The impact scope represents the radial distance when the temperature difference is close to zero. Fig. 8 indicates that the impact scope with full heat insulation is larger than that without heat insulation, which can eventually reach 13 m after 120 days. And it is evident that the temperature gradient with full thermal insulation has a much more dramatic change near the wellbore, which is beneficial to heat extraction. Furthermore, as a contrast, the production time has little influence on the impact scope without thermal insulation, meaning its thermal process achieves stable faster. But this would lead to a lower heat output. Then the temperature difference distribution with full insulation along axial direction is analyzed in Fig. 9. It is clear the impact scope and temperature gradient increase with well depth. The reason is that a larger temperature difference between working fluid and the deep reservoir would enhance heat transfer.
ACCEPTED MANUSCRIPT
Fig. 8. Comparison of temperature difference distribution with full thermal insulation and without thermal insulation at 30d, 60d, 90d, 120d
Fig. 9. Temperature difference distribution with full thermal insulation at 1550 m, 2550 m, 3550 m, 4550 m To determine the impact scope and further analyze the sustainability of CBHE, the temperature difference distribution before the next space heating in buildings is analyzed in Fig. 10, just 8 months in length. Apparently, the temperature difference distribution in axial direction has the
ACCEPTED MANUSCRIPT similar trend changing with that in Fig. 9. Compared to the temperature difference at 4550 m after the first space heating, the temperature difference at 4550 m in Fig. 10 even decreases by 92.63 °C. This means that the temperature field in formation almost recovers its original condition, verifying the sustainability of CBHE. In addition, note that the impact scope has an increase of 14 m, which shows that the determination of the well spacing needs to consider the effect of recovery period. Hence, further study is carried out in Fig. 11. It is clear the impact scope increases with recovery period, which ranging from 21 m to 35 m. Furthermore, the decrease trend of temperature difference at 4550 m tends to be slow, and it only has maximum decline of 2.57 °C from 8th month to 20th month. This also indicates that there is no need to employ a longer recovery period for the next indoor heating.
Fig. 10. Temperature difference distribution at 1550 m, 2550 m, 3550 m, 4550 m after 8 months of recovery
ACCEPTED MANUSCRIPT
Fig. 11. Temperature difference distribution at 4550 m under different recovery periods 5. Analysis of heat extraction performance To achieve the best heat extraction performance, it is highly necessary to study the effects of the key factors, including inlet temperature, mass flow, insulation section length of central tube, thermal conductivity of cement and reservoir conditions, etc. And the formula for thermal power is as follows:
Qout 0.001q out Cout Tout in Cin Tin
（20）
5.1 The effect of flow rate Based on the ground surface condition, the inlet temperature is set to 25 °C. Fig. 12 illustrates that the heat output decreases drastically in a particularly short time because of the sharp change of temperature near the wellbore. Additionally, it is difficult to compensate the energy loss around the borehole in time. Hence, the decreasing rate of outlet temperature would slow down gradually, then remains relatively stable. And the outlet temperature declines with flow rates, which leads to a smaller difference between the outlet temperature and inlet temperature in Eq. 20. But due to the increase of the fluid mass participating in heat transfer in unit time, the thermal power would have an upward trend. In view of the heat supply requirement, the outlet temperature of CBHE needs to be kept above 60 °C to maintain the appropriate indoor temperature [34]. This means that the flow rate should be less than 25 m3/h in Fig. 13, and its outlet temperature is approximately 61.81 °C at 120d. Furthermore, the change flow rate has a dramatic effect on the pressure drop, which is a quite important parameter during production. Hence, the pressure drop under different flow rates is investigated in Fig. 14. And we employ the assumption of a smooth tube for the calculation, which has been proved its good precision for a U-tube collector [35]. It is clear the pressure drop has a noticeable increase with flow rate, which can range from 3.729 MPa to 22.021 MPa. Moreover, the growth rate of pressure drop increases as rises the flow rate, while that of thermal power slows down
ACCEPTED MANUSCRIPT gradually. Therefore, in terms of pressure drop and thermal power, 25 m3/h would be best for the CBHE performance. Due to the pressure drop at 25 m3/h satisfies the production requirements, so the flow rate of 25 m³/h is employed for the following analysis.
Fig. 12. Outlet temperature and thermal power with production time under different flow rates
Fig. 13. Outlet temperature and thermal power with flow rate at 30d, 60d, 90d, 120d
ACCEPTED MANUSCRIPT
Fig. 14. Thermal power and pressure drop with flow rate 5.2 The effect of inlet temperature Fig. 15. illustrates that the outlet temperature increases with inlet temperature, but the thermal power declines linearly, and its curve slope is approximately -13.55 kW/°C at 120d. This indicates that the inlet temperature has significant effect on the heat extraction performance of CBHE. And the lower temperature difference between the reservoir and working fluid would weaken the heat transfer. Accordingly, if the outlet temperature has meet the heating requirements, there is no need to raise the inlet temperature. Thus, the 25 °C serve as the inlet temperature in following study based on the surface condition.
Fig. 15. Outlet temperature and thermal power with inlet temperature at 30d, 60d, 90d, 120d 5.3 The effect of insulation section length Fig. 16 shows that the growth rate of outlet temperature and thermal power increases when the
ACCEPTED MANUSCRIPT insulation length is less than 3500 m, then this trend becomes slow gradually and remains basically unchanged after it exceeds 4500 m, verifying the feasibility of partial insulation. And its thermal power has only a drop of 1.12% compared to that with full insulation at 120d. The reason is that the temperature difference between annulus and central pipe declines gradually with well depth according to the temperature distribution of working fluid in Fig. 7, which can weaken the thermal process near the bottom. Additionally, the difference of outlet temperature and thermal power tends to be larger with insulation length. By comparison with the production curve in Fig. 12, it can be found that a shorter insulation length would make the production stable in less time and further lead to a drop of total energy output.
Fig. 16. Outlet temperature and thermal power with insulation section length at 30d, 60d, 90d,120d 5.4 The effect of cement sheath Because of the lower thermal conductivity of cement compared to that of formation, so it may hinder the heat transfer between the reservoir and working fluid. This phenomenon is very similar to the skin effect in petroleum engineering, which has significant influence on the oil production [36, 37]. Obviously, the outlet temperature and thermal power have a significant change when the thermal conductivity is less than 1.0 W/(m·K) in Fig. 17, then this trend slows down gradually, verifying the existence of skin effect in geothermal development. And the thermal power at 2.8 W/(m·K) (corresponding to the formation) can increase by 12.22% compared to that at 0.7 W/(m·K) after 120 days. Note that if the well cementing operations is not carried out, the formation water (0.569~0.687 W/(m·K) [38]) would fill its space, which can result in a worse performance. Apparently, the heat extraction region mainly lies in the deep formation, so there is no need to employ high thermal conductivity cement in whole length. Fig. 18 displays that the outlet temperature difference and thermal power difference with high thermal conductivity length. This is conducted by subtracting the outlet temperature and thermal power at a certain high thermal conductivity length from that with high thermal conductivity cement in whole length, respectively. The thermal conductivity of cement sheath is set to 2.8 W/(m·K). It should be emphasized that 0m on the x axis means the depth of well bottom, while 3760 m on the x axis represents the depth of
ACCEPTED MANUSCRIPT intermediate casing. It can be seen the trends of reduction of outlet temperature difference and thermal power difference tend to be slow with high thermal conductivity length. And there is only an outlet temperature difference of 0.31 °C at 3760 m after 120 days, which proves the feasibility of partial cementing with high thermal conductivity cement.
Fig. 17. Outlet temperature and thermal power with thermal conductivity of cement at 30d, 60d, 90d,120d
Fig. 18. Outlet temperature difference and thermal power difference with high thermal conductivity length at 30d, 60d, 90d,120d 5.5 The effect of geological conditions The geological conditions are particularly important basis for the location optimization of a geothermal well, so the effects of the thermal conductivity of formation and geothermal gradient is analyzed. Fig. 19 shows that the outlet temperature and thermal power increase with the thermal
ACCEPTED MANUSCRIPT conductivity of formation, and this trend slows down gradually. Furthermore, the thermal power at 4.0 W/(m·K) ((corresponding to granite[39]) is approximately 1290.24 kW compared to the 555.29 kW at 1.0 W/(m·K) (corresponding to shallow soil[39]) after 120 days, and its percentage of increase reaches 132.35%. This also indicates that the CBHE system is more suitable for the development of deep geothermal resources. Fig. 20 illustrates that the outlet temperature and thermal power increase linearly with geothermal gradient, and the curve slope can even reach 13.06 °C/(°C/100m) and 379.19 kW/(°C/100m) at 120d, respectively. Consequently, the thermal conductivity of formation and geothermal gradient both have remarkable impact on the heat extraction performance.
Fig. 19. Outlet temperature and thermal power with thermal conductivity of formation at 30d, 60d, 90d,120d
Fig. 20. Outlet temperature and thermal power with geothermal gradient at 30d, 60d, 90d,120d 5.6 The effect of well depth
ACCEPTED MANUSCRIPT For the cost assessment of a geothermal well, the well depth is a very important factor [40], so the thermal characteristics under different well depth are investigated to provide guidance for the well depth optimization. Fig. 21 depicts that the growth rate of outlet and thermal power increase with well depth because of the heat transfer enhancement between the reservoir and working fluid. And if we employ a well depth of 2000 m, the outlet temperature difference can reach 30.88 °C compared to the outlet temperature at 5000 m (corresponding the design depth). Hence, the CBHE system has great superiority in exploiting the deep geothermal resources. Furthermore, it would take longer to achieve stable as rises the well depth according to the production curve in Fig. 12, which can further lead to the increase of total energy output.
Fig. 21 Outlet temperature and thermal power with well depth at 30d, 60d, 90d,120d 6. Conclusion In this paper, an unsteady-state heat transfer model is proposed for a deep CBHE. The finite difference method is employed to solve the heat transfer model, and the model is validated with experiment data. Subsequently, the temperature field in formation is investigated. And the influences of key parameters on heat extraction performance of CBHE are analyzed. The main conclusions are as follows: (1) The temperature field analysis shows that there is a significant temperature gradient around the borehole, far more than the geothermal gradient, indicating the heat transfer mainly occurs in the radial direction. The impact scope in the reservoir is relatively limited, so the CBHE has a great potential in developing geothermal resources. Furthermore, the temperature field in formation before the next heating period would get close to its initial condition, which illustrates that the heat production of CBHE system is sustainable. And the impact scope increases with recovery period, which is an important reference for well spacing optimization. (2) The outlet temperature and thermal power have a remarkable decrease at the initial stage, but then remains relatively stable. Moreover, a higher flow rate leads to the decline of outlet temperature, but the thermal power would increase. Further study indicates that in terms of the pressure drop and thermal power, the 25 m3/s is the optimal flow rate which could obtain the best CBHE performance. (3) As the inlet temperature increases, the outlet temperature declines and the thermal power increases. This depicts that a larger temperature difference between the reservoir and working fluid
ACCEPTED MANUSCRIPT could enhance heat transfer. Therefore, a lower inlet temperature is recommended if the outlet temperature has meet the required temperature. (4) The temperature of working fluid in central tube only has a little reduction, verifying its well heat insulation effect. Additionally, the insulation section length of central tube only has a significant effect on heat extraction performance within a certain range. And the changes of outlet temperature and thermal power are minor when this length exceeds 4500 m, which indicates that there is an optimal insulation length. (5) A lower thermal conductivity of cement sheath could weaken the heat transfer from the reservoir to working fluid, verifying the existence of the skin effect in geothermal development. Consequently, it is necessary to conduct the research of high thermal conductivity cement. And to reduces its cost, further analysis depicts that there is an optimum length with high thermal conductivity cement because the thermal process mainly lies in the deep reservoir. (6) The increase of the thermal conductivity of formation and geothermal gradient can greatly enhance the heat exchange, so it is particularly important to determine the appropriate geological conditions for CBHE system. Additionally, the well depth has a noticeable influence on the thermal process, and even the growth rate of outlet temperature and thermal power increase with well depth. The findings can offer guidance for the location optimization and well depth optimization of CBHE system. Acknowledgements The authors would like to acknowledge the National Key Research and Development Program of China (Grant No. 2016YFE0124600). Besides, support from the Program of Introducing Talents of Discipline to Chinese Universities (111 Plan) (Grant NO. B17045) is appreciated. Nomenclature AS = annular space a = geothermal gradient, °C/m CBHE = coaxial borehole heat exchanger CP = central pipe Cpm = heat capacity of working fluid, J/(kg·°C) Cin = heat capacity of inlet fluid, J/(kg·°C) Cout = heat capacity of outlet fluid, J/(kg·°C) Cw = heat capacity of casing, J/(kg·°C) Cc = heat capacity of cement sheath, J/(kg·°C) Cs = heat capacity of formation, J/(kg·°C) Deq = hydraulic diameter, m f = Darcy friction factor h1 = convective heat transfer coefficient for inner surface of central pipe, W/(m2·°C) h2 = convective heat transfer coefficient for outer surface of central pipe, W/(m2·°C) h3 = convective heat transfer coefficient for casing wall, W/(m2·°C) Nu = Nusselt number Pr = Prandtl number Qout = thermal power, kW q = flow rate of working fluid, m³/h R = thermal resistance, (m·°C)/W Re = Reynolds number
ACCEPTED MANUSCRIPT r1 r2 r3 r4 r5 r6 r7 Tf1 Tf2 Tin Tout Tw Tc ρm ρin ρout ρw ρc ρs λf λw1 λw2 λw3 λw λc λwc λs λcs Δt Δz i j
= = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
inner radius of central pipe, m inner radius of insulation layer, m outer radius of insulation layer, m outer radius of central pipe, m inner radius of casing, m outer radius of casing, m outer radius of cement sheath, m temperature of working fluid in central pipe, °C temperature of working fluid in annulus space, °C temperature of inlet fluid, °C temperature of outlet fluid, °C temperature of casing, °C temperature of cement sheath, °C density of working fluid, kg/m3 density of inlet fluid, kg/m3 density of outlet fluid, kg/m3 density of casing, kg/m3 density of cement sheath, kg/m3 density of formation, kg/m3 thermal conductivity of working fluid, W/(m·°C) thermal conductivity of inner layer of central pipe, W/(m·°C) thermal conductivity of insulation layer of central pipe, W/(m·°C) thermal conductivity of outer layer of central pipe, W/(m·°C) thermal conductivity of casing, W/(m·°C) thermal conductivity of cement sheath, W/(m·°C) equivalent thermal conductivity of casing and cement, W/(m·°C) thermal conductivity of formation, W/(m·°C) equivalent thermal conductivity of cement and formation, W/(m·°C) time step, s length of axial grid, m number of grid cells in radial direction number of grid cells in axial direction
References [1] J.C. Reboredo, Renewable energy contribution to the energy supply: Is there convergence across countries?, Renewable & Sustainable Energy Reviews 45 (2015) 290-295. [2] M. Sun, P.P. Zhang, T.H. Shan, C.C. Fang, X.F. Wang, L.X. Tian, Research on the evolution model of an energy supply– demand network, Physica A Statistical Mechanics & Its Applications 391(19) (2012) 4506-4516. [3] M.E. Bildirici, S.M. Gökmenoğlu, Environmental pollution, hydropower energy consumption and economic growth: Evidence from G7 countries, Renewable & Sustainable Energy Reviews 75 (2017) 68-85. [4] Z. Lyu, X. Song, G. Li, X. Hu, Y. Shi, Z. Xu, Numerical Analysis of Characteristics of a Single U-tube Downhole Heat Exchanger in the Borehole for Geothermal Wells, Energy 125 (2017) 186-196. [5] R.A. Beier, Transient heat transfer in a U-tube borehole heat exchanger, Applied Thermal Engineering 62(1) (2014) 256-266. [6] H. Holmberg, J. Acuña, E. Næss, O.K. Sønju, Thermal evaluation of coaxial deep borehole heat exchangers, Renewable Energy
ACCEPTED MANUSCRIPT 97 (2016) 65-76. [7] A. Zarrella, M. Scarpa, M.D. Carli, Short time-step performances of coaxial and double U-tube borehole heat exchangers: Modeling and measurements, Hvac & R Research 17(6) (2011) 959-976. [8] R.N. Horne, Design considerations of a down-hole coaxial geothermal heat exchanger, Trans. - Geotherm. Resour. Counc.; (United States) 4 (1980). [9] P.J. Yekoladio, T. Bello-Ochende, J.P. Meyer, Design and optimization of a downhole coaxial heat exchanger for an enhanced geothermal system (EGS), Renewable Energy 55(4) (2013) 128-137. [10] R.A. Beier, J. Acuña, P. Mogensen, B. Palm, Transient heat transfer in a coaxial borehole heat exchanger, Geothermics 51(7) (2014) 470-482. [11] D. Gordon, T. Bolisetti, S.K. Ting, S. Reitsma, Short-term fluid temperature variations in either a coaxial or U-tube borehole heat exchanger, Geothermics 67 (2017) 29-39. [12] D. Gordon, T. Bolisetti, S.K. Ting, S. Reitsma, A physical and semi-analytical comparison between coaxial BHE designs considering various piping materials, Energy 141 (2018). [13] K. Morita, T. Yamaguchi, H. Karasawa, H. Hayamizu, Development and evaluation of temperature simulation code for geothermal wells, J. Min. Metall. Inst. Jpn 100 (1984) 1045-1051. [14] K. Morita, W.S. Bollmeier, H. Mizogami, Analysis of the results from the Downhole Coaxial Heat Exchanger (DCHE) experiment in Hawaii,
(1992).
[15] M.S. Saadi, R. Gomri, Investigation of dynamic heat transfer process through coaxial heat exchangers in the ground, International Journal of Hydrogen Energy 42(28) (2017). [16] E. Zanchini, S. Lazzari, A. Priarone, Improving the thermal performance of coaxial borehole heat exchangers, Energy 35(2) (2010) 657-666. [17] E. Zanchini, S. Lazzari, A. Priarone, Effects of flow direction and thermal short-circuiting on the performance of small coaxial ground heat exchangers, Renewable Energy 35(6) (2010) 1255-1265. [18] X. Song, Y. Shi, G. Li, R. Yang, Z. Xu, R. Zheng, G. Wang, Z. Lyu, Heat Extraction Performance Simulation for Various Configurations of a Downhole Heat Exchanger Geothermal System, Energy 141 (2017). [19] Y. Shi, X. Song, G. Li, R. Li, Y. Zhang, G. Wang, R. Zheng, Z. Lyu, Numerical Investigation on Heat Extraction Performance of a Downhole Heat Exchanger Geothermal System, Applied Thermal Engineering
(2018).
[20] Y. Fujimitsu, K. Fukuoka, S. Ehara, H. Takeshita, F. Abe, Evaluation of subsurface thermal environmental change caused by a ground-coupled heat pump system, Current Applied Physics 10(2) (2010) S113-S116. [21] K. Bär, W. Rühaak, B. Welsch, D. Schulte, S. Homuth, I. Sass, Seasonal High Temperature Heat Storage with Medium Deep Borehole Heat Exchangers ☆, Energy Procedia 76 (2015) 351-360. [22] G. Cui, S. Ren, L. Zhang, J. Ezekiel, C. Enechukwu, Y. Wang, R. Zhang, Geothermal exploitation from hot dry rocks via recycling heat transmission fluid in a horizontal well, Energy 128 (2017) 366-377. [23] Morita, I I u AN EXPERIMENT TO PROVE THE CONCEPT OF THE DOWNHOLE COAXIAL HEAT EXCHANGER (DCHE) IN HAWAII, Geothermal Resources Council Trans 16 (1992) 9-16. [24] T. Kohl, M. Salton, L. Rybach, Data analysis of the deep borehole heat exchanger plant weissbad (Switzerland), Proc World Geothermal Congress
(2000).
[25] T. Kohl, R. Brenni, W. Eugster, System performance of a deep borehole heat exchanger, Geothermics 31(6) (2002) 687-708. [26] J. Acuña, B. Palm, Distributed thermal response tests on pipe-in-pipe borehole heat exchangers, Applied Energy 109(6) (2013) 312-320. [27] J. Zhu, K. Hu, X. Lu, X. Huang, K. Liu, X. Wu, A review of geothermal energy resources, development, and applications in China: Current status and prospects, Energy 93 (2015) 466-483. [28] J.W. Lund, T.L. Boyd, Direct utilization of geothermal energy 2015 worldwide review, Geothermics 60 (2016) 66-93. [29] J. Hou, M. Cao, P. Liu, Development and utilization of geothermal energy in China: Current practices and future strategies,
ACCEPTED MANUSCRIPT Renewable Energy 125 (2018). [30] A.D. Chiasson, G.G. Culver, D. Favata, S. Keiffer, DESIGN, INSTALLATION, AND MONITORING OF A NEW DOWNHOLE HEAT EXCHANGER, IEEE International Conference on Granular Computing, 2005. [31] V. Gnielinski, New equations for heat and mass transfer in the turbulent flow in pipes and channels, Nasa Sti/recon Technical Report A 75(2) (1975) 8-16. [32] G.K. Filonenko, Hydraulic Resistance in Pipes, Teploenergetika 4(4) (1954) 15-21. [33] Y. Shi, X. Song, G. Li, Z. Shen, X. Hu, Z. Lyu, R. Zheng, G. Wang, Numerical Analysis of the Heat Production Performance of a Closed Loop Geothermal System, Renewable Energy 120 (2017). [34] T. Nussbaumer, S. Thalmann, Status report on district heating systems in IEA countries, Swiss Federal Office of Energy, and Verenum, Zürich
(2014).
[35] J. Acuña, Improvements of U-pipe Borehole Heat Exchangers, Annan Naturresursteknik (2010). [36] A.F.V. Everdingen, The Skin Effect and Its Influence on the Productive Capacity of a Well, Journal of Petroleum Technology 5(6) (1953) 171-176. [37] W. Hurst, J.D. Clark, E.B. Brauer, The skin effect in producing wells, Journal of the Society of Petroleum Evaluation Engineers (1968). [38] E.O. Holzbecher, Modeling Density-Driven Flow in Porous Media, Springer1998. [39] G. Hellström, Ground heat storage: Thermal analyses of duct storage systems, Physics & Astronomy
(1991).
[40] LUKAWSKI M. Z, ANDERSON B. J, AUGUSTINE Chad, CAPUANO L. E, BECKERS K. F, LIVESAY Bill, TESTER J. W, Cost Analysis of Oil, Gas, and Geothermal Well Drilling, Journal of Petroleum Science & Engineering 118 (2014) 1-14.
ACCEPTED MANUSCRIPT
Highlights of “Numerical analysis of heat extraction performance of a deep coaxial borehole heat exchanger geothermal system” Xianzhi Song1*, Gaosheng Wang1, Yu Shi1, Ruixia Li2, Zhengming Xu1, Rui Zheng1, Yu Wang1, Jiacheng Li1 (1. State Key Laboratory of Petroleum Resources and Prospecting, China University of Petroleum, Beijing, Beijing 102249, China; 2. Sinopec Star Petroleum Co. Ltd, Beijing 100083, China)
An unsteady-state heat transfer numerical model for CBHE is presented. The simulation values show particularly good agreement with field data. The temperature field dynamic response of CBHE is analyzed comprehensively. Influences of key factors on the heat production performance of CBHE are studied.