Paper - Stanford School of Earth Sciences

PROCEEDINGS, Thirty-Ninth Workshop on Geothermal Reservoir Engineering
Stanford University, Stanford, California, February 24-26, 2014
SGP-TR-202
Temperature Behavior of Geothermal Wells During Production, Injection and Shut-in
Operations
Kaan Kutun, Omer Inanc Tureyen and Abdurrahman Satman
ITU Maden Fakultesi, Petrol ve Dogal Gaz Muh. Bol., Maslak 34469 Istanbul, TURKEY
[email protected], [email protected], [email protected]
Keywords: Temperature profile, static conditions, dynamic conditions, numerical modeling
ABSTRACT
Static and dynamic measurements of pressure and temperature along geothermal wells are commonplace practices for
characterization purposes in the geothermal industry. Such profiles provide good insight regarding the deliverability of the well,
location of the upper and lower boundaries of the reservoir / reservoirs and etc. The static temperature profiles, for example, can be
used for modeling the natural state of the system through history matching. As these profiles are taken, it is important that the actual
static and dynamic conditions are reached. In other words, to take a dynamic temperature profile, one would have to wait until the
temperature profile stabilizes in the well after production has started and vice versa for the static profile. Hence it becomes crucial
to know beforehand how long the stabilization time takes. In this study a numerical model is developed to study the parameters
effecting the stabilization time of static and dynamic conditions in the wells for a single phase water system. The model is based on
mass and energy balance equations and couples the reservoir with the well and takes into account the heat losses to the
surroundings of the well. The model is validated using various analytical models in the literature. A synthetic application is
provided to identify key parameters that effect temperature distributions. In the application we model temperature behavior along
the well for the transition from static to dynamic and from dynamic to static.
1. INTRODUCTION
Static and dynamic temperature profiles taken along geothermal wells provide very useful information regarding the fluid and
petrophysical properties of the reservoir. For example, dynamic temperature profiles could help identify different zones of water
entry (at different temperatures) into the well. Dynamic temperature profiles taken at different times could give how the bottomhole
temperature changes with time. Or they could provide insight into how much heat is lost to the surroundings of the well as the fluid
is moving in the well. Static profiles on the other hand could be used to identify the top and bottom points of the reservoir itself. In
fact multiple entry points may be determined from the static profiles.
When interpreting either the static or dynamic temperature profiles along the wells, two points are very crucial to consider. The first
is that a mathematical model which properly describes the physics of the phenomenon must exists so that inference of various fluid
and/or petrophysical properties can be made. The second is that during the measurements, a long enough time must be considered
for either the static or dynamic conditions to be established. Failing to do so could result in a mischaracterization of the system. The
necessary time required to reach static or dynamic conditions can be assessed through the use of an appropriate mathematical
model.
Many authors have developed models that describe temperature behavior in geothermal wells for various cases. Perhaps the earliest
and the most cited work is given by Ramey (1962). In Ramey’s work, an approximate solution to the wellbore heat transmission
problem for the injection of hot or cold fluids is provided. This solution assumes steady state flow in the wellbore and heat transfer
into the surrounding formations is treated to be unsteady radial conduction. Wooley (1980) presents a numerical model which
models the heat conduction effects in the horizontal and vertical directions as well as modeling the convective nature of the flow.
Another numerical model is given by Farouq Ali (1981) which presents a comprehensive numerical model that can deal with
different well operation conditions. This model is capable of handling steam/water mixtures and is based on solving the momentum,
mass balance and energy balance equations in the wellbore. Durrant and Thambynayagam (1986) give a straightforward iterative
procedure for the upward and downward movement of a steam/water mixture. Wu and Pruess (1990) have presented an analytical
solution that considers heat losses to an arbitrary number of layers with different properties. The simplifying assumptions of Ramey
(1962) were not included in this study. Hagoort (2004) assess Ramey’s (1962) method for the calculation of temperatures in
injection and production wells. This study states that although Ramey’s approximation and a rigorous solution are in good
agreement at later times, Ramey’s model overestimates the temperatures at early transient periods. Hasan and Kabir (2010) give a
model for two phase flow using the drift-flux approach. Livescu et. al. (2010) give a semi analytical approach where the extension
of isothermal wellbore-flow models to non-isothermal cases are presented.
2. THE MATHEMATICAL MODEL
The mathematical model developed in this study is based on solving the mass balance and energy balance equations simultaneously
for a given set of grid blocks. The structure of the model is very similar to that given by Tureyen and Akyapi (2011). Lets consider
any grid i given in Figure 1. We assume that the grid is composed of water and rock components. Grid i can make an arbitrary
number of connections with any other grid in the system. The total number of connected grids is termed Nc. Energy and mass
transfer is allowed for between grid i and the connected grids represented by the index j. The mass balance applied to grid i is given
in Eq. 1.
1
Kutun et al.
d   w
V b ,i
i
N c ,i

dt
   p
i , jl
jl

 pi  
w ,i , jl
z
jl
  w
 zi
p ,i
 w inj , i  0
(1)
l 1
Here, Vb denotes the bulk volume,  the density,  the porosity, t time, α the transmissibility,  the pressure gradient, p the pressure
and w the mass rate. The subscripts w represents water and b represents bulk. The first term of Eq.1 represents the mass
accumulation. The second term represents the sum of the mass transfer between grid i and the connected grids and takes into
account the flow component due to gravity. The third term gives the mass production rate and the fourth term gives the mass
injection rate. Here it is important to note that the mass injection rate is performed at a specified temperature Tinj. The mass transfer
between the grids is based on the pressure difference between the grids and the transmissibility. The transmissibility between the
grid block i and any neighboring grid block jl can be defined as given in Eq. 2.
 kA  

 d   i , jl
 i , j  
l
(2)
Here k is the permeability, A is the cross-sectional area, d is the distance between the grid points of grids i and jl and µ is the
viscosity. Equation 2 gives the transmissibility term for two neighboring grid blocks in Cartesian coordinates. Different expressions
for the transmissibility may be obtained for different coordinate systems such as a radial coordinate system or a spherical
coordinate system.
Grid : j1
 i, j
1
 i, j
Wp,i
Production
Ti
Winj,i
Injection
Tinj,i
Grid : j2
2
 i, j
Grid : i
Water + rock
 i, j
Grid : j3
3
l
Grid : jl
Volume : Vbi
Porosity : i
Temperature : Ti
Pressure : pi
 i, j
N c i 1
Grid : jNci-1
 i, j
N ci
Grid : jNci
Figure 1: Illustration of any grid i and the neighboring grids.
The energy balance equation used in the model is given in Eq. 3.
V b ,i
d
dt
1     m C m T
N c ,i


l 1
h  i , j
l
 p
jl
  w u w i   wh


inj
  wh
w
i
N c ,i
 p i   i, j z j  zi
l
w
l
    T
i , jl
jl
 Ti

(3)
l 1
Here C is the heat capacity, u is the internal energy, h is the enthalpy and  is the heat conduction transmissibility. The subscript m
represents the matrix. The first term in Eq. 3 is the accumulation of energy in grid block i. As it is clear, the accumulation of energy
takes place in both components of the grid block; both in the rock and in the water. The second term represents the energy
contribution due to injection, the third term represents the energy contribution due to production, the fourth term represents the
convective heat transport from and to the neighboring grid blocks and the final term represents the conductive heat transfer to and
from the neighboring grid blocks. Here it is important to note that an upwinding scheme is applied to the convective heat transfer
between grid blocks. The subscript  denotes the direction of upwinding and is defined in Eq. 4.
2
Kutun et al.
 i if p
  

j l if p
jl
jl
 pi  
 pi  
w ,i , jl
z
z
w ,i , jl
jl
jl

 z 
 zi
(4)
i
The heat conduction transmissibility for Cartesian coordinates is defined as given in Eq. 5. Just as in evaluating the flow
transmissibility term given in Eq. 2, the heat conduction transmissibility term given in Eq. 5, can also be written for different
coordinate systems.
 A 
i , j  

l
 d  i , jl
(5)
Here  represents the thermal conductivity.
In evaluating the thermal conductivity and the permeability for the interface between grid blocks i and jl in Eq’s 2 and 5, harmonic
averages are used. For the other parameters evaluated at the interface (such as the density or viscosity) arithmetic averaging is used.
The mass and energy balance equations are treated in a fully implicit manner causing them to become highly non-linear. Hence a
Newton Raphson procedure is used to solve Eq’s 1 and 3 simultaneously. For constructing the Jacobian matrix in the Newton
Raphson procedure, numerical derivatives are used. A forward difference scheme is used to handle the derivatives with respect to
time.
The above mathematical method has been used previously by Palabiyik et. al. (2013) and later again by Palabiyik (2013) for a
detailed sensitivity analysis of factors effecting pressure and temperature behaviors in geothermal reservoirs. In this study we
extend this model to model the wellbore temperature behavior. The schematics of the grid blocks used in this model are given in
Figure 2.
Strata grid blocks
=0
Well grid blocks
=1
α=0
=computed
Specified Production or injection
α=0
=computed
Reservoir grid blocks
Figure 2: Schematics of the grid blocks used in the study.
As can be seen in Figure 2, the well is also discretized in the z direction using the same discretization as the reservoir and the
overburden strata. Afterwards, the mass and energy balance equations are solved in the wellbore grids just as we do so in the
reservoir and strata grids. The only difference is that a much higher permeability is used for the wellbore grids. Throughout this
study in all our runs, the permeability of the wellbore grid blocks are taken to be 110-5 m2. This kind of an approach where we use
Darcy type of a flow for the wellbore allows us to couple the reservoir flow with flow in the wellbore without solving the
momentum balance equation in the wellbore. To eliminate the rock component from the wellbore grid blocks, the porosity of the
wellbore is taken to be 1. The thick lines separating the reservoir grids from the strata grids represent the reservoir upper boundary.
At this boundary the flow transmissibility is set to zero for avoiding flow into the strata grids. The heat conduction transmissibility
however is set to a computed value to allow for heat transfer from the reservoir into the strata by way of conduction. The same
approach is taken for the boundary between the well grids and the strata grids. This way heat loss to the formations by way of
3
Kutun et al.
conduction can be modelled while preventing fluid movement into the strata grids. Mass production or injection rates are specified
only at the very top grid of the wellbore. No other mass rate is specified for any other grid block.
3. VERIFICATION OF THE MODEL
In this section we provide a verification of the model using two of the common analytical approaches for modeling the temperature
behavior of wells; the approach given by Ramey (1962) and the approach provided by Hagoort (2004). Table 1 lists the properties
used in the example for the verification case. Figure 3 gives the grid system used in the verification case. The initial distribution of
temperature has been obtained using a temperature gradient of 0.118 C /m. Initially the same temperature gradient is assumed to
have established in the well. Only one grid layer in the z direction is used to represent the reservoir. The rest of the grids in the z
direction are used for the strata. For the verification, we consider only a production scenario.
Table 1: Well and reservoir properties for the verification example.
Well radius, m
Reservoir radius, m
Reservoir permeability, m2
Reservoir porosity, fraction
Rock compressibility, bar-1
Rock thermal expansion coefficient, C-1
Rock density, kg/m3
Rock specific heat capacity, J/kg-C
Thermal conductivity of water, J/m-s-C
Thermal conductivity of rock, J/m-s-C
Depth of reservoir, m
Height of reservoir, m
Initial pressure (@ 1050 m), bar
Initial temperature (@ 1050 m), C
Number of grid blocks in r direction
Number of grid blocks in z direction
Production rate, kg/s
0.1
1000
110-12
0.2
2.910-5
0
1952
1000
0.67
2.92
1000
100
117.85
144.09
11
15
4.39
Well
Strata
Reservoir
Figure 3: The structure of the grids used in the verification example.
The results of the verification example are given in Figure 4. In Figure 4 we provide a comparison of well head temperature
between our model, Ramey’s (1962) model and the model proposed by Hagoort (2004). To perform the comparison, the solutions
presented in Figure 5 of Hagoort (2004) are used. The Ramey number and the Graetz number for these solutions are one. Hence, to
do the comparison, the values presented in Table 1 give the same Ramey and Graetz numbers. The solutions have been digitized
and compared with the model presented in this work in Figure 4. As can be seen in Figure 4, the rigorous solution given by Hagoort
and the model presented in this study match rather well both for early times and for late times. Matches for early times cannot be
4
Kutun et al.
obtained with the Ramey model since the Ramey model tends to overestimate the temperatures for early times as stated by Hagoort
(2004).
140
130
120
Ramey
Hagoort
Model
110
Temperature, C
100
90
80
70
60
50
40
30
20
10
Time, Days
Figure 4: Comparison of various models.
4. SYNTHETIC APPLICATIONS
In thıs section we provide various synthetic applications of the model specifically geared towards analyzing the parameters
affecting the stabilization time of the temperature. We consider a production case and a shut in case and point out the most
important parameters that effect the stabilization time. In all the synthetic applications, presented in this section, the parameters
given in Table 1 are used unless otherwise stated.
Production Case
In this subsection, we analyze with more detail the synthetic application for the verification example given in the previous section.
We first provide how the temperature profile inside the well evolves with time. Figure 5 gives the temperature profiles for various
time slices of 110-6, 0.01, 0.1, 1, 10, 100 and 1000 days. It is important to note that, most of the change occurs during the first 10
days. Then the changes in profile with time become very small. At 110-6 days, there is practically no change in the temperature
profile. The distribution is almost identical to the initial temperature distribution in the well. As time progresses, the linear behavior
of the profile is distorted and the profile starts to shift. Two main mechanisms of heat transfer play a role in the changing well
temperature profile. The first is the convective heat transfer from the bottom of the well to the top via production of the fluid. The
second is the conductive heat loss to the surroundings of the well. In the final profile at 1000 days, the wellhead temperature
reaches to about 133 C different from the initial reservoir temperature which was at around 144 C. This difference between the
wellhead temperature and the bottomhole temperature is because of the heat losses to the surroundings of the well.
Next we consider the effect of the production rate on the stabilization time of the temperature. However, it is first important to
discuss briefly the physics of the problem. From a purely mathematical point of view, it is never possible to reach a stabilized
temperature in the well simply because of the conductive heat losses to the surroundings of the well. However, from a practical
point of view, at late times during the production, the temperature does not change much with time as we have seen in Figure 5. In
order to analyze better the change of temperature with time, we look at the derivative of the wellhead temperature with respect to
time. Figure 6 gives the behavior of this derivative with time. As can be seen from Figure 6, the derivative initially displays a
constant behavior. Then, it decreases more or less linearly. This is indicative of two different behaviors. The constant derivative
portion of the curve reflects the time period where convective flow is dominating the temperature change. The linearly decreasing
portion of the data reflects the period where the conductive heat losses now dominate the temperature change.
In order to compare the stabilization times, an arbitrary derivative value of 0.1 C/s for the wellhead temperature is chosen as the
point where stabilization is said to have occurred. In other words, it is assumed that the temperature has stabilized if the derivative
decreases below this cut off during production.
5
Kutun et al.
0
20
40
Temperature, C
60
80
100
120
140
160
0
100
200
300
Depth, m
400
500
600
700
800
900
1000
t=0.000001 Days
t=0.01 Days
t=0.1 Days
t=1 Days
t=10 Days
t=100 Days
t=1000 Days
1100
Figure 5: Temperature profiles of a producing well.
10000
1000
dT/dt, C/D
100
10
1
0.1
0.01
0.001
Time, Days
Figure 6: The behavior of the temperature derivative with time.
We first compare the stabilization times for various mass flow rates. The results are given in Figure 7. As it is clear, the
stabilization times decrease hyperbolically with increasing mass flow rate. There are two reasons for this. The first reason is
because at lower mass flow rates it takes a longer time for the hotter fluid to travel to the wellhead. The second reason is that at
6
Kutun et al.
lower mass flow rates, as the fluid is moving towards the wellhead, more heat is lost to the surroundings of the well, causing the
fluid to arrive at the wellhead with lower temperatures compared to what would have been with higher mass flow rates.
Next we look at the effect of the well radius on the stabilization time. Figure 8 illustrates the results. An increase in the well radius
causes an increase in the stabilization time. This is because, at a constant mass flow rate an increase in the radius results in a
decrease in the velocity of the fluid. Hence it takes a longer time for the hotter fluid to arrive at the wellhead.
50
45
Stabilization Time, Days
40
35
30
25
20
15
10
5
0
2
4
6
8
10
12
14
Mass Flow Rate, kg/s
16
18
20
0.16
0.18
0.2
Figure 7: Wellhead temperature stabilization times for various mass flow rates.
32
30
Stabilization Time, Days
28
26
24
22
20
18
16
14
0.02
0.04
0.06
0.08
0.1 0.12 0.14
Well Radius, m
Figure 8: Wellhead temperature stabilization times for various well radii.
7
Kutun et al.
Shut In Case
In this section we only consider the effects of production time on the stabilization time of temperature. For this purpose we perform
production for various durations of time. Then the well is shut in and the stabilization times are observed. The same cut-off
derivative value of 0.1 is used for the shut-in period. Figure 9 gives the results. According to Figure 9 as the production increases
the stabilization time during shut in increases as well. After a certain point, the stabilization time becomes constant. As mentioned
earlier, during production the surroundings of the well are heated due to the hot fluid flowing in the well. Larger production times
causes the heated region around the well to become wider. Hence more heat is stored in the surroundings of the well. This leads to
longer stabilization times during shut in since it takes a longer time for the surroundings of the well cool with more heat stored
during production.
105
100
Stabilization Time, Days
95
90
85
80
75
70
65
60
55
0
100
200
300 400 500 600 700
Production Time, Days
800
900 1000
Figure 9: Wellhead temperature stabilization times for various production durations.
5. CONCLUSIONS
The following conclusions have been obtained from this study:






A model capable of modeling temperature profiles in the well is developed. The developed model is also coupled to the
reservoir allowing for realistic profiles to be computed.
The developed model has been verified by two of the common methods used in the literature; the method of Ramey
(1962) and the rigorous solution given by Hagoort (2004).
The focus of this study are the stabilization times of wellhead temperature for various flow rates, well radii and producing
times for shut in.
It is found that with increasing flow rates, the stabilization times for temperature during production are decreased.
An increase in the well radius causes an increase in the wellhead temperature stabilization times during production.
Increasing the production time causes an increase in the stabilization time of the wellhead temperature for shut in.
REFERENCES
Durrant A. J. and Thambynayagam, R. K. M.: Wellbore Heat Transmission and Pressure Drop for Steam/Water Injection and
Geothermal Production: A Simple Solution Technique, SPE Reservoir Engineering, 1, (1986), 148-162.
Farouq Ali, S. M.: A Comprehensive Wellbore Steam/Water Flow Model for Steam Injection and Geothermal Applications, SPE
Journal, 21, (1981), 527-534.
Hagoort, J.: Ramey’s Wellbore Heat Transmission Revisited, SPE Journal, 9, (2004), 465-474.
Hasan, A.R. and Kabir, C. S.: Modeling Two-Phase Fluid and Heat Flows in Geothermal Wells, Journal of Petroleum Science and
Engineering, 71, (2010), 77-86.
Livescu, S., Durlofsky, L. J. and Aziz, K.: A Semianalytical Thermal Multiphase Wellbore-Flow Model for Use in Reservoir
Simulation, SPE Journal, 15, (2010), 794-804.
Palabiyik, Y: A Study on Pressure and Temperature Behaviors of Geothermal Wells in Single-Phase Liquid Reservoirs, PhD
Dissertation, Istanbul Technical University Graduate School of Science Engineering and Technology, Istanbul, Turkey (2013).
8
Kutun et al.
Palabiyik, Y., Tureyen, O. I., Onur, M. and Paker Deniz, M: Pressure and Temperature Behaviors of Single-Phase Water
Geothermal Reservoirs Under Various Production/Injection Schemes, Proceedings, 38th Workshop on Geothermal Reservoir
Engineering, Stanford University, Stanford, CA, USA (2013).
Ramey, H.J.: Wellbore Heat Transmission, Journal of Petroleum Technology, 14, (1962), 427-435.
Tureyen, O.I., and Akyapi, E.: A Generalized Non-Isothermal Tank Model for Liquid Dominated Geothermal Reservoirs,
Geothermics, 40, (2011), 50-57.
Wu, Y.S. and Pruess, K.: An Analytical Solution for Wellbore Heat Transmission in Layered Formations, SPE Reservoir
Engineering, 5, (1990), 531-538.
9

Download Report