Investigation of thermal model effects on geothermal and oil reservoir well testing behavior during water and steam injection

. Global warming and reducing fossil fuel resources have increased the interest in using renewable resources such as geothermal energy. In this paper, in the ﬁ rst step, heat transfer equations have been presented for reservoir during water (steam) injection by considering heat loss to adjacent formations. According to radius of thermal front, the reservoir is partitioned into two regions with different ﬂ uid physical properties. The heat transfer model is coupled with a ﬂ uid ﬂ ow model which is used to calculate the reservoir pressure or ﬂ uid ﬂ ow rates. Then by calculating outer radius of heated region and using radial composite reservoir model, the ﬂ uid ﬂ ow equations in porous media are solved. Using pressure derivative plot in regions with different thermal conductivity coef ﬁ cients, a type curve plot is presented. The reservoir and adjacent formation thermal conductivity coef ﬁ cients can be calculated by matching the observed pressure data on the thermal composite type curve. Additionally, the interference test in composite geothermal reservoir is discussed. In the composite reservoir model, parameters such as diffusivity coef ﬁ cient, conductivity ratio and the distance to the radial discontinuity are considered. New type curves are provided to introduce the effect of diffusivity/conductivity contrast ratios on temperature behavior. Improving interpretations, and performing fast computations and fast sensitivity analysis are the bene ﬁ ts of the presented solutions.


Introduction
With a decline in conventional oil resources and rising demand, extraction of heavy oil resources has become more prominent (Speight, 2016). Heavy oil is characterized with high viscosities and low API gravities less than 20 (Attanasi and Meyer, 2010). Heavy oil is widely distributed in many countries (Butler, 1991). The most important way to extract heavy oil is to reduce the viscosity of oil and hence improve oil mobility (Ali, 2003). Thermal recovery techniques, such as Cyclic Steam Injection (CSI), Steam Flooding (SF), in situ combustion and Steam Assisted Gravity Drainage (SAGD), are the most widely used techniques to develop heavy oil reservoirs (Upreti et al., 2007). Cyclic Steam Injection (CSI) process has been the primary thermal technique for the extraction of heavy oil reservoirs in the past decades (Al Adasani and Bai, 2011;Kovscek, 2012). A cyclic steam-injection method has three stages. The first stage of this method is the slug injection of steam into the reservoir. In the second stage of cyclic steam-injection method, the well is shut in for several days to render a uniform heat distribution. In the third stage, oil is produced from the same injection well. This method is repeated as long as it has economic justification (Speight, 2016). After cyclic steaminjection method, the reservoir is divided into two finite zones with different temperature and viscosity, heated zone and cold zone (Green and Willhite, 1998).
Several models have been developed to predict cyclic steam injection process in vertical wells. Fontanilla and Aziz (1982) wellbore model can predict bottom hole pressure. Marx and Langenheim (1959) described a method for predicting thermal investigation radius and cumulative heated area for hot fluid injection at a constant rate. Boberg and Lantz (1966) developed the first analytical model of cyclic steam injection for a vertical well whose radius of the steam zone is based on Marx and Langenheim model. Bentsen and Donohoe (1969) developed a transient flow model to calculate oil rate.
Generally, well tests from thermal recovery processes are analyzed using a radial, two region composite model (Ambastha, 1989;Ambastha and Ramey Jr., 1989). The study performed by Theis (1935) was the first article which presented an analytical solution for semi-infinite aquifer with constant rate. Perina and Lee (2006) developed a general semi analytical model for constant pumping in a finite aquifer with a partially penetration well and skin. Wang and Yeh (2008) gave a review on the relationship between the transient and steady-state solution for constant rate and pressure tests in finite and infinite aquifers. Wang et al. (2012) presented an analytical model for predestining the drawdown pressure due to a constant rate from a well with a skin zone in finite aquifers. Su et al. (2015) developed three-region composite transient pressure analysis model for CO 2 flooding by considering the skin and well bore storage effect of well. Huang et al. (2015) developed a new solution for modeling steady-state radial and vertical flows in a composite aquifer. Lin et al. (2016) developed an analytical model for predicting transient drawdown pressure induced by pumping at a radial composite finite aquifer with Robin-type boundary condition. Cheng et al. (2001) simulated the heat extraction by water injection from a fracture embedded in geothermal reservoir. Boyadjiev et al. (2005) extended the Lauwerier formulation for describing temperature of porous media saturated with oil. Stopa and Wojnarowski (2006) presented an analytical solution to obtain the front velocity of injected water in a geothermal reservoir. Ghassemi et al. (2008) investigated poro-thermoelastic effects of water injection in an enhanced fracture geothermal system by considering fluid flow through fracture in a hot rock. Shaik et al. (2011) investigated the importance of heat transfer between the matrix block and fluid on economic hot water production from a naturally fractured geothermal reservoir. Shaw-Yang and Hund-Der (2008), Li et al. (2010) and Yeh et al. (2012) presented semi-analytical models by using the Laplace transform method in order to model the temperature distribution in an aquifer thermal extraction system for a finite aquifer with different rock properties. Ascencio et al. (2014) developed an analytical model for cold water injection into a fractured reservoir which presents solutions to compute the distance that separates the hydrodynamic and thermal fronts within the reservoir at a given time under transient and steady state conditions. Ganguly and Kumar (2014) presented an analytical model to predict the transient temperature distribution during cold water injection in a heterogeneous geothermal reservoir. Abbasi et al. (2018) attempted to develop a generalized analytical solution to Advection-Conduction Equation (ACE) for the phenomenon of heat transfer through porous media considering transient boundary conditions. They additionally evaluated heat transfer shape factor during both pseudo-steady state and transient periods, and further validated their analytical solution with numerical finite difference and element schemes.
The work performed by App (2013) put forth the effect of hydraulic fractures on flowing bottom-hole temperatures under 1000-5000 psi drawdown conditions. It was illustrated that due to the flow regime change from radial for a non-hydraulically fractured completion to linear for a hydraulically fractured completion, the Joule-Thomson expansion thermal effect is reduced. App (2016) scrutinized the flow geometry effect on sandface temperature for high drawdown circumstances where the impact of Joule-Thomson thermal is tangible. The temperature characteristics for both linear and radial regimes were examined when the flow was either transient or steady-state. Via examining the thermal-energy equation from a dimensionless analysis standpoint, the authors presented a mechanistic evaluation of the terms introduced.
Al Saedi et al. (2018) put forward novel analytical fluidtemperature models with disparate well bottom hole boundary conditions by formulating second-order differential equation that can capture more factual downhole conditions. The incentive to develop his methodology (which could handle the casing diameter) was the lack of published field thermal data in the literature. Chevarunotai et al. (2018) introduced an analytical approach to assess temperature distribution in a single-phase oil reservoir assuming constant rate production, along with an application of the presented method in single-phase gas reservoirs. In addition, App's approach was used in the presented model to incorporate heat transfer from a reservoir to both underburden and overburden formations. Hashmi et al. (2015) highlighted two objectives in their work: 1to assess rate from Multipoint-Temperature Sensor (MTS) data available for downhole telemetry with a rigorous transient analytical model, and to show the possible use of the computed rates to perform transient analysis for both the pre-and post-cleanup periods. Their methods entailed transient-temperature modeling with and without superposition effects.
In this paper, initially by using the heat transfer equations, the temperature distribution plot versus the distance from the wellbore is calculated. Then, by using viscosity correlations the viscosity profile plot versus temperature is obtained. Using pressure derivative plot in different thermal conductivity coefficients, a type curve plot is presented. The reservoir and adjacent formation thermal conductivity coefficients can be calculated based on well pressure data and the mentioned type curve. Also, in this paper, the interference test in composite geothermal reservoir is discussed. All important parameters are considered for the composite reservoir model such as conductivity ratio, diffusivity coefficient and the distance to the radial discontinuity. In the end, type curves are presented.

Mathematical model
In this section, we will describe the analytical model of temperature distribution when cold water is injected into geothermal reservoirs (steam in oil reservoirs) and then, based on the creation of two regions with different physical properties, we will describe fluid flow equations for a composite reservoir. In the following, an analytical solution of the heat transfer model will be expressed in a composite reservoir.

Reservoir heating model
The heat transfer equation for fluid flow through a heavy oil reservoir (or a geothermal reservoir) is considered in this study, where the reservoir is assumed to be confined by impermeable rock layers. A schematic representation of the system with an injection well is given in Figure 1.
Assumptions that are used in developing the model are as follows: 1. The heavy oil reservoir (geothermal reservoir) is of infinite radius and with homogeneous properties, which is bounded above and below by impermeable formation. All the rock and fluid properties in the reservoir are assumed constant. (Li et al., 2010;Stopa and Wojnarowski, 2006). 2. The rock layers thickness are assumed to be uniform and homogeneous. 3. The heat transfer from the reservoir to the rock media is assumed to be one-dimensional and horizontal thermal conduction is neglected. 4. The initial temperature prior to the injection is assumed to be constant. Geothermal gradient varies between 0.6 and 1.6°F/100 ft and hence assuming total of e.g. 300 ft pay thickness, only 3°F temperature difference exists from the top most to lower most layer. 5. Due to large surface area exist in a porous media through which the heat transfer can take place and also the large heat capacity of fluids (hydrocarbon/ water) in place, thermal equilibrium established with no sensible lateral temperature gradient.
By applying the principle of energy conservation together with the above assumptions, the governing equations for temperature distributions in the aquifer and adjacent formation can be formulated as (Gringarten and Sauty, 1975;Li et al., 2010;Stopa and Wojnarowski, 2006): where, For the upper/lower formation, the heat conduction equation can be written as, The analytical solution of equation (1) was evaluated by Li et al. (2010), who presented a model for predicting the transient temperature distribution. The reservoir temperature in Laplace space becomes: where,

Reservoir composite model
The mathematical model for composite reservoir with its schematic diagram shown in Figure 2, is subject to the following fundamental assumptions: 1. The formation is homogeneous, isopachous, and isotropic. 2. The reservoir is filled with slightly compressible fluid. 3. Flow is isothermal and follows Darcy's law. 4. The gravity and capillary force are negligible.
The dimensionless diffusivity equation for two different regions: where subscripts 1 and 2 denote the heated zone and unheated zone, respectively, r D is the dimensionless radial distance, R D is the outer radius of the heated zone, t is the production time, and g is Initial condition: Inner boundary conditionsconstant rate: If inner boundary condition is constant pressure, then: Outer boundary condition: In order to satisfy continuity for the pressure and flux at the interface between the heated zone and unheated zone, the following equations should be considered, respectively: The dimensionless variables are defined as follows: Dimensionless pressure for constant rate: Dimensionless pressure for constant pressure: By using the Laplace transforms method, the solutions of equations (4a) and (4b) subject to equations (5a)-(5g) were evaluated by Satman et al. (1980), who presented a model for a radial composite under constant rate in a finite-radius well. Also, Turki et al. (1989) presented a model for a radial composite under pressure.

Temperature profile in composite geothermal reservoir
The problem studied in this part is the heat transfer in a composite geothermal reservoir. The region r w r < r a is of one property and region r a r < 1 of another. The governing equations are: For region-1: For region-2: where a ¼ k qc . The initial temperature is: Inner boundary condition: Outer boundary condition: Inter boundary conditions: In order to reduce the complexity, it is preferred to solve the equation in dimensionless format. For this purpose, the dimensionless parameters are defined as follows: ð9fÞ Substituting equations (9a)-(9g) into equations (7a) and (7b) and (8a)-(8e) yields: The initial is: Inner boundary condition: Outer boundary condition: Inter boundary conditions: Using the same procedure as Turki et al. (1989), the solution to equations (10a) and (10b) with boundary conditions can be expressed in Laplace domain as follows: where

Results and discussion
At first the reservoir temperature at each point of the reservoir is calculated using heat transfer model which can be applied for cold water injection into geothermal reservoirs or steam injection into heavy oil reservoirs. The graphs of pressure versus time are obtained from equations (8a) to (8e) and the effect of different reservoir parameters on them will be considered later. Also, based on equations (11a) and (11b), the interpretations of temperature response at the observation well during water injection are presented for composite geothermal reservoirs.

Steam injection in heavy oil reservoir
Steam injection is one of the production methods from heavy oil reservoirs which contain high viscosity oil. In this type of reservoirs, steam is cyclically injected for a few months then oil is produced from the same well. During the steam injection due to heat loss to underlying and overlying rock layers, heat front warms up the reservoir and reduces oil viscosity. After steam injection, the oil reservoir is divided into two zones with different viscosity values; the near wellbore region that heat front of steam has reached there and the second zone that still remains at initial reservoir temperature. By using the presented information in Table 1 and equation (3), reservoir temperature can be calculated at each time and location.
To validate the derived analytical solution, a radial model was numerically solved. As shown in Figure 3, cell coordinates were specified in terms of r-h-z values and log spacing was used in the radial direction to better simulate and capture the wellbore response. In order to verify the one-dimensional flow, the permeability is set to zero in the z and h directions. The rock and fluid properties used to simulate Cyclic Steam Injection (CSI) model are given in Table 1 and the relative permeabilities and capillary pressure were taken from SPE Case 1.
The reservoir oil temperature versus location is shown in Figure 4. As seen in this figure, the near wellbore area has the maximum temperature and due to heat front advancement and heat loss through adjacent formation layers, the temperature decreases and finally reaches to initial reservoir temperature. The result of the Cyclic Steam Injection (CSI) analytical model show a good agreement with the simulation result.
The oil viscosity versus location is shown in Figure 5. The Beggs and Robinson (Beggs and Robinson, 1975) oil viscosity correlation is used to calculate viscosity respect to temperature.
The average temperature in Region-1 can be calculated as: where r a is radius of thermal front. The average temperature of region 1 versus time is shown in Figure 6. As shown in the figure, with increasing the duration of steam injection in the reservoir, the average temperature of region 1 increases. According to Figure 7, after cyclic steam injection into oil reservoir, the thermal front of the steam has progressed to a distance from the wellbore which causes oil viscosity alteration and the reservoir is partitioned into two regions with different mobility ratios. The result of the Cyclic Steam Injection (CSI) analytical model presented in this work and those based on simulation tally up.
In Figure 8, the flow rate response of the oil reservoir after Cyclic Steam Injection (CSI) process is presented. The upper curve in Figure 8 demonstrates the response of wellbore flow rate. The response of early time portion is infinite acting, and because the outer region mobility is lower, the discontinuity acts as a no-flow boundary and the rate declines exponentially. At late times, the flow rate is controlled by characteristics of the outer region and the response is infinite acting. The rate response at the radial discontinuity versus time is shown in the lower curve in Figure 8. Initially the discontinuity flow rate is zero, and increases with time as fluid flows from the outer region to the inner region. At late time, the rate at the heat front is practically the same as wellbore rate and the flow rate of inner region is constant. k h parameter represents the relative thermal conductivity coefficient of adjacent formation to the reservoir. As shown in Figure 9, by increasing the k parameter, the thermal front radius will be reduced. This means that increasing the thermal conductivity of adjacent formation causes higher rate of heat transfer to adjacent formation layers. As a result, heat front progresses into smaller radius around the wellbore and the altered viscosity region diminishes.
Pressure derivative plot versus different values of k h is shown in Figure 10. Pressure derivative plot for production well after cyclic steam injection shown in Figure 10 that is similar to well test plots for composite reservoirs. Also, the dual porosity/dual permeability models exhibit similar derivative fingerprint when the early fracture flow is very short or it is masked with the wellbore phenomena. This may cause reservoir characterization misinterpretation. As seen in Figure 10 and in the interpretation of Figure 9, by increasing the value of k h parameter, radius of the heated region is reduced and the composite reservoir pattern in pressure derivative plot happens earlier. On the other hand, according to expressed pattern plot in Figure 10, the thermal conductivity coefficient of reservoir and adjacent formation can be calculated by matching well pressure/production data to pressure derivative plot.

Water injection in composite geothermal reservoir
For specific values of dimensionless radial discontinuity, dimensionless distance, diffusivity ratio, and conductivity ratio the numerical inversions are prepared. The dimensionless temperature response at the observation well in Region-2 is investigated under the following conditions: (1) varying dimensionless radial discontinuity and dimensionless distance to the observation well, (2) varying diffusivity ratio and conductivity ratio of the composite system. These two cases are considered independently,  and the effects of various dimensionless parameters on the temperature response at the observation well are studied. The dimensionless temperature response at the observation well, T 2D , is computed as a function of dimensionless time, t 2D , for various values of dimensionless radial discontinuity, R D , and for constant values of dimensionless distance, r D , diffusivity ratio, n, and conductivity ratio, j. In Figure 11 T 2D is shown as a function of t D for different values of R D . Figure 12 shows the dimensionless temperature response versus (t D ) for different values of j and constant value of (R D ). As shown in Figure 12, by increasing the conductivity      ratio (decreasing the thermal conductivity of Region-2 or increasing the thermal conductivity of Region-1) the Region-2 acts like an isolator for the Region-1 and the effect of injected water in composite reservoir is arrived later to the observation well and eventually the value of temperature reduction will be higher.

Conclusion
In the present communication, heat transfer equations were presented for reservoir during steam injection process by considering the effects of heat loss to adjust formations. According to the radius of thermal front, the reservoir was divided into two distinct thermal regions. Afterward the fluid flow equations in porous media was solved by calculating radius of the heated region and applying the well test equations of a composite reservoir. The thermal conductivity coefficient of the reservoir and adjust formation were calculated based on the well pressure data and presented type curve which is pressure derivative plot in different thermal conductivity coefficients. A novel solution for analyzing the temperature response at the observation well in a composite geothermal reservoir was presented here. Using the obtained solution, the temperature response at the observation well can be expressed as a function of dimensionless radial discontinuity and conductivity ratio. Moreover, the type-curve matching performed with the obtained data from the observation well provides practical information regarding composite reservoir properties. The log-log approach for determining the properties of Region-2, needs more analysis time compared to the typecurve matching. Improving the interpretations, fast sensitivity analysis and fast computations are the benefits of the derived solutions.