Static Formation Temperature Prediction Based on Bottom Hole Temperature

Abstract

Static formation temperature (SFT) is required to determine the thermophysical properties and production parameters in geothermal and oil reservoirs. However, it is not easy to determine SFT by both experimental and physical methods. In this paper, a mathematical approach to predicting SFT, based on a new model describing the relationship between bottom hole temperature (BHT) and shut-in time, has been proposed. The unknown coefficients of the model were derived from the least squares fit by the particle swarm optimization (PSO) algorithm. Additionally, the ability to predict SFT using a few BHT data points (such as the first three, four, or five points of a data set) was evaluated. The accuracy of the proposed method to predict SFT was confirmed by a deviation percentage less than ±4% and a high regression coefficient R² (>0.98). The proposed method could be used as a practical tool to predict SFT in both geothermal and oil wells.

1. Introduction

Deep drilling is necessary for the exploitation of deep geothermal reservoirs [1]. In this case, borehole drilling is a complicated process in which a constant thermal anomaly (in addition to the circulating drilling mud) affects the static formation temperature (SFT) around the borehole [2]. Determining SFT at any depth demands a lot of time to measure the bottom-hole temperature (BHT) and shut-in time [3]. Measuring BHT can be costly due to the usage of sophisticated logging equipment and the necessity to temporarily stop the wellbore drilling [4].

Optimal estimation of SFT is required for several applications, including the determination of geothermal heat flow, analysis of well logs, estimation of geothermal potential, evaluation of in situ thermophysical formation properties [5], and the determination of hydrocarbon properties in petroleum systems [6,7,8].

Estimation of SFT is usually achieved by analytical and numerical simulation methods. Most of the analytical methods are based on the constant linear and cylindrical heat source models. Analytical methods most commonly used to estimate SFT include the Horner-plot method (HM), or the line-source method [9]; the Kutasov-Eppelbaum method (KEM), or the generalized Horner method [10]; the Manetti method (MM), or the cylindrical source with a conductive heat flow method [11]; the Hansan and Kabir method (HK), or the cylindrical heat source with a conductive-convective heat flow method [12]; the Bernnand method (BM), or the radial source with a conductive heat flow method [13]; the spherical and redial heat flow method (SRM) proposed by Ascencio et al. [14]; and the Leblanc method (LM), or the cylindrical source with a conductive heat flow method [15].

These methods determine SFT by using BHT and shut-in time data as input, and the linear or nonlinear regression models as solutions [16]. Nevertheless, large errors are likely encountered in the prediction of SFT. In this case such errors may arise from various sources, including unrealistic models proposed to describe the drilling process, heat transfer models based on simple assumptions, measurement errors in the BHT data, and total uncertainties in SFT estimation [17].

Numerical simulation is another method to estimate SFT, which can also be applied to determine geothermal gradients and describe the thermal history [18]. For example, Garcia et al. [19] developed the numerical simulator TEMLOPI for estimating the transient temperature distribution in a wellbore and the surrounding rock formation. Application to well Az-29 from the Los Azufres Mexican geothermal field shows satisfactory results. This simulator could be used by drilling engineers to determine the optimal design of cement slurries and their setting times during well construction.

However, both analytical and numerical simulation methods have some limitations, such as the excessive amount of some other data that is needed in addition to BHT and shut-in time, (e.g., thermophysical and transport properties of the wellbore, drilling and cementing materials, and formation and rock materials, whose data is rarely available and, therefore, limits the usage of these methods) and the accurate circulation time that is usually unknown or difficult to determine under drilling conditions.

Against this background, a reliable and practical tool to estimate SFT is still required in geothermal and petroleum industries. In this paper, a mathematical function has been proposed to correlate BHT and shut-in time. The coefficients of the function were obtained from the particle swarm optimization (PSO) algorithm based on the least squares fit target.

PSO is a stochastic, population-based optimization method that was introduced by Kennedy and Eberhart [20]. It belongs to the family of swarm intelligence computational techniques and is inspired by social interaction in human beings and animals (especially bird flocks and fish schools). PSO optimizes a problem by having a population of candidate solutions, dubbed particles here, and moving these particles around in the search-space according to simple mathematical formulae over the particle’s position and velocity. Each particle’s movement is influenced by its local best known position, but is also guided toward the best known positions in the search-space, which are updated as better positions are found by other particles. The PSO algorithm has been used in many numerical solution problems and shows its wide applicability [21,22,23].

PSO has some advantages in solving optimization problems, for example, it requires few parameters to be tuned by users, it is highly accurate, it is less affected by initial solutions compared with other algorithms, it has fast convergence, and it is easy to code due to the simple underlying concepts and demands no requirement for preconditions, such as continuity or differentiability of the objective functions [24].

2. Methodology

A function correlating BHT and shut-in time was derived to fit the BHT data and the estimated SFT. The coefficients of this function can be obtained from least squares fit method using the particle swarm optimization (PSO) algorithm. Additionally, other methods were also introduced for comparison purposes with the new method. Statistical tests were applied to evaluate the validity of the predicting methods.

2.1. Method Development

2.1.1. Function Derivation

The Horner method for determining the static formation temperature is widely used in the oil and gas industry [9]. This analytical method is based on the assumption that thermal effect of drilling is a constant linear heat source. The approximate solution is given by:

$$\text{BHT}\left(\mathrm{t}\right)={\mathrm{T}}_{HM}-\left({\mathrm{b}}_{HM}\right)·\text{ log}\left\{\left({\mathrm{t}}_C+\mathrm{t}\right)/\mathrm{t}\right\}$$

(1)

where ${\mathrm{T}}_{HM}$ is the static formation temperature, $\text{ log}\left\{\left({\mathrm{t}}_C+\mathrm{t}\right)/\mathrm{t}\right\}$ is the dimensionless Horner time (DHT), and ${\mathrm{t}}_C$ and $\mathrm{t}$ are the circulation time before shut-in and the time elapsed since the circulation stopped, respectively. One of the problems in Equation (1) is that of being inconsistent with the boundary conditions. BHT approaches the static formation temperature, ${\mathrm{T}}_{HM}$ when t approaches to infinity. However, BHT cannot be obtained when t approaches or is equal to zero. This problem may also decrease the fitting quality of the model. We propose the following mathematical model in order to solve this problem in Equation (1). The modified model is expressed as:

$$\text{BHT}\left(\mathrm{t}\right)=\mathrm{a}-\mathrm{b}·\text{ log}\left\{\frac{1}{\left(1+\mathrm{c}·\mathrm{t}\right)}+1\right\}.$$

(2)

where a, b, and c are constants, a is actually equal to ${\mathrm{T}}_{HM}$ , which can be estimated using Equation (2) when shut-in time tends to infinity:

$$\widehat{\text{SFT}}=\underset{\mathrm{t}\to \infty }{\lim }\text{BHT}=\mathrm{a}$$

(3)

Equation (2) meets all of the boundary (time) conditions:

As t approaches infinity, maximum BHT is obtained:

$${\text{BHT}}_{max}={\mathrm{T}}_{HM}=\mathrm{a}$$

(4)

When t approaches to zero, the minimum BHT is obtained:

$${\text{BHT}}_{min}=\mathrm{a}-\mathrm{b}·\text{ ln}2$$

(5)

One can observe that the problem in Equation (1) has now been solved by Equation (2). When the maximum and minimum values of BHT have been determined, the shape of the BHT-time curve will only depend on the value of c. The curve of the BHT-time function is illustrated in Figure 1. The equation can characterize the BHT-time function in a large scope, as shown in Figure 1.

2.1.2. Solutions to the Three Parameters in the New Function

In order to obtain the best fit of the equation, the least squares fit target is applied. The difference between the proposed function and the measured value of temperature is denoted by Q, which is obtained using the following equation:

$$\mathrm{Q}={\sum }_{i=1}^n{\left({\text{BHT}}_i-{\text{BHT}}_{ci}\right)}^2$$

(6)

where ${\text{BHT}}_i$ and ${\text{BHT}}_{\text{ci}}$ are measured and calculated BHT using Equation (2).

The least squares fit method requires minimization of Q. Here we used the particle swarm optimization (PSO) algorithm to obtain the best fitting coefficients of a, b, and c that yields the minimum value of X.

The position and velocity of the ith particle are respectively denoted by $X_i$ and $V_i$, which are the best position and velocity attained by each particle and therefore referred to as its personal best. For the ith particle, the position vector of its personal best is denoted by $P_i$ and its objective value is denoted by $P_{best}$. The best position attained by the swarm is referred to as either the swarm best, global best, or leader. It is denoted by $P_g$ and its objective is denoted by $g_{best}$. At each iteration t, the positions and velocities of all particles are updated by Equations (7) and (8) which are also called the update equation [24]:

$$V_i\left(t+1\right)=\omega V_i\left(t\right)+C_1{r}_1\left(P_i-X_i\right)+C_2{r}_2\left(P_g-X_i\right)\left(i=1,2,\dots ,N_P\right)$$

(7)

$$X_i\left(t+1\right)=X_i\left(t\right)+V_i\left(t+1\right)\left(i=1,2,\dots ,N_P\right)$$

(8)

where $\omega$ is inertia weight; $C_1$ is the cognitive acceleration coefficients, which prompts the attraction of the particle towards its own personal best; $C_2$ is the social acceleration coefficients, which prompts the attraction of the particle towards the swarm best; and $r_1$ and $r_2$ are two random numbers in [0,1].

A schematic of the solving process of the proposed method is depicted in Figure 2.

2.2. Existing Methods

2.2.1. Selection of the Analytical Methods

A few previous analytical methods have been discussed in this paper for comparison with the newly-proposed method in this paper. The three most commonly used analytical methods were introduced in

Table 1. A summary of the analytical methods used for comparison with the method developed in this paper.

Method	Physical Model	Equation	Sources
Horner (HM)	Constant Linear heat source	$\text{BHT}\left(\mathrm{t}\right)={\mathrm{T}}_{HM}-{\mathrm{b}}_{HM}ln\left(\frac{t+t_c}{t}\right)$	Dowdle and Cobb (1975)
Manetti (MM)	Conductive cylindrical heat source	$\text{BHT}\left(\mathrm{t}\right)={\mathrm{T}}_{MM}-{\mathrm{b}}_{MM}ln\left(\frac{t}{t-t_c}\right)$	Manetti (1973)
Ascencio (SRM)	Spherical-radial heat flow	$\text{BHT}\left(\mathrm{t}\right)={\mathrm{T}}_{SRM}-{\mathrm{b}}_{SRM}ln\left(\frac{1}{\sqrt{t}}\right)$	Ascencio et al. (1994)

. The equations, analytical models and sources of these methods are listed in Figure 1. Basically, two regression models were applied to obtain the coefficients in each analytical method.

2.2.2. Selection of the Regression Models

Two regression methods were used as regression models to obtain the coefficients in the three analytical methods: the ordinary linear regression (OLR) model and quadratic regression (QR) model.

The general equation of the OLD model is: $\mathrm{y}=\mathrm{a}+\text{bx}$ , where a and b are the intercept and the slope of the fitted straight line [16], while that for the QR model is given by: $\mathrm{y}=\mathrm{a}+\text{bx}+c{x}^2$ , where a, b, and c represent the polynomial coefficients [16].

2.3. Statistical Evaluation

2.3.1. Deviation Percentages (DEV%)

In order to test the estimation accuracy of the estimated SFT, the deviation percentages (DEV%) between the estimated and reference SFT values were used:

$$\text{DEV}\%=[\frac{\left(\widehat{SFT}-SFT\right)}{SFT}]\times 100$$

(9)

where $\widehat{SFT}$ is the value of SFT estimate, and SFT is the true SFT value reported in the literature.

It should be noted that this evaluation test was not applied from datasets without reference SFT (Data 7 and Data 8). The closer the DEV% to zero, the better the estimation of SFT.

2.3.2. Regression Coefficient (R²)

The regression coefficient (R²) is applied for testing the fitting ability of each method. The value can be obtained from:

$$R^2=1-\frac{{\sum }^{\text{}}{\left(y_i-{\widehat{y}}_i\right)}^2}{{\sum }^{\text{}}{\left(y_i-\overline{y}\right)}^2}.$$

(10)

where $y_i$ is the measured value of BHT, ${\widehat{y}}_i$ is the estimated BHT, and $\overline{y}$ is the mean of the measured BHT values. The fitting is more satisfactory if the value of R² is closer to 1.

2.3.3. Residual Sum of Squares (RSS)

The fitting result of each method can also be evaluated through estimation of the normalized residual sum of squares (RSS), which is calculated by the following equation:

$$\text{RSS}=\frac{{\sum }_{i=1}^n{\left(y_i-{\widehat{y}}_i\right)}^2}{n}.$$

(11)

where $y_i$ and ${\widehat{y}}_i$ are same as in Equation (10), n is the total numbers of elements in a BHT dataset. The smaller the value of RSS in a dataset, the better the fitting result. The goodness of fit can be evaluated by DEV% and R². However, the comparison between different methods should be based on the whole accuracy and fitting ability of each method. To evaluate each method comprehensively, the following synthetic statistical parameters were also used.

2.3.4. Theil Inequality Coefficient (TIC)

TIC is applied to test the estimation accuracy comprehensively when only very few BHT data are available. For each dataset and each method used, the TIC value can be obtained from:

$$\text{TIC}=\frac{\sqrt{{\sum }_{i=3}^5{\left({\widehat{SFT}}_i-SFT\right)}^2}}{\sqrt{{\sum }_{i=3}^5{{\widehat{SFT}}_i}^2}+SFT}$$

(12)

where ${\widehat{SFT}}_i$ is the SFT estimate using the first $\mathrm{i}$ of the dataset.

The value of TIC is in [0,1], and could be used to compare between different datasets. For each method, the closer the TIC to zero, the more accurate the method is.

2.4. Data Sources

Eight thermal recovery datasets were collected from the published literature (see

Table 2. Summary of the bottom-hole temperature (BHT) datasets used in this paper.

Data	Type	n	tc (h)	Sources	Data Name in This Paper
SHBE	Synthetic data	8	5	Shen and Beck (1986)	Data 1
CLAH	Synthetic data	15	5	Cao et al. (1988)	Data 2
CJON	Synthetic data	12	0.2	Cooper and Jones (1959)	Data 3
KJ-21	Geothermal field data	6	2.5	Steingrimsson and Gudmundsson (1989)	Data 4
SG	Geothermal field data	12	3	Schoeppel and Gilarranz (1966)	Data 5
MOU	Synthetic data	3	10	Mou (2013)	Data 6
DA-XIN	Geothermal field data	40	5	Da-Xin (1986)	Data 7
UASM	Petroleum field data	14	10	Kutasov (1999)	Data 8

) for the accuracy and application tests. They include:

(1)

Four synthetic datasets selected from the literature; and

(2)

Four datasets logged in some boreholes from long logging work in geothermal and petroleum fields.

These datasets are summarized in Table 2 [25,26,27,28,29,30].

It should be pointed out that the datasets 1–6 all have reference SFT values, which can be very useful to evaluate the application of the proposed method and conduct comparisons between different methods.

3. Validation and Discussion

3.1. Accuracy of the Static Formation Temperature (SFT)

For each dataset, coefficients a, b, and c were solved by the PSO algorithm based on the least squares fit target. After applying Equation (2), the BHT data were recalculated at each shut-in time of the borehole and entitled calculated BHT. The DEV% between SFT estimate and reference SFT were obtained by Equation (9), as described in the “Methodology” section. A comparison between data from the literature and data estimated by our proposed method is depicted in Figure 3, in which yellow circles represent measured BHT values (in this case reference BHT), black dots represent calculated BHT, red full lines represent reference SFT, and the black dashed lines represent values of SFT estimate.

Using the synthetic set of Data 1, the SFT value estimated by the proposed method was 80.99 °C, which is in agreement with the true SFT, 80 °C. The error in this case is 1.24%. The SFT estimated in Data 2 is about 4 °C higher than the true value (120 °C), with an acceptable percentage of 3.52%. For those data sets from well-logging (Data 4 and 5), the deviation percentages fall in the range of −2% to 1%. Data 6 was recorded at three shut-in times, with an SFT value of 105.296 °C estimated accurately by the proposed method, and a deviation percentage of 0.28% with the reference SFT value of 105 °C.

Figure 4 shows the comparison between estimated SFT values using the proposed method and true SFT values for six datasets using the true SFT values ranging from 20.25 to 240 °C, during the last shut-in time from 1.5 h to 600 h. The number of the dataset varies from three to 15. The DEV% of each dataset is in [−2%, 4%], a satisfactory range for a practical predicting tool.

3.2. Fitting Results

The regression coefficients R² determined by Equation (10) and the normalized RSS by Equation (11) were calculated for all datasets using the proposed method. The calculation results are presented in

Table 3. Results of each dataset using the proposed method.

NO.	Coefficients of the Proposed Method			STF	SFT	DEV%	R2	RSS
	a	b	c	(Reference)	(Estimated)
Data 1	81.000	45.960	0.172	80.000	81.000	1.249	1.000	0.025
Data 2	124.232	63.459	0.243	120.000	124.232	3.527	0.999	0.089
Data 3	20.327	30.865	28.299	20.250	20.327	0.380	0.987	0.073
Data 4	236.084	280.837	0.039	240.000	236.084	−1.632	0.999	2.243
Data 5	100.132	21.111	0.694	101.111	100.132	−0.968	0.999	0.026
Data 6	105.296	26.481	0.067	105.000	105.296	0.282	1.000	0.000
Data 7	132.459	133.588	0.289	N/A	132.459	N/A	0.979	11.053
Data 8	147.600	12.809	0.071	N/A	147.600	N/A	0.979	0.085

, and the BHT values calculated by the proposed method are also shown in Figure 3. For the datasets 1 and 2, the values of R² are both greater than 99.9%, which indicates that the proposed method can fit the BHT data accurately. As for those BHT datasets obtained from well-logging, the proposed method also showed a good fitting ability because the R² for each dataset is greater than 98%. Based on the good matching results, the proposed method seems to provide acceptable correlations between BHT and shut-in time.

In addition, comparisons of estimates and reference SFT values using various methods was carried out, the results of which are shown in

Table 4. Comparison of estimates and reference Static formation temperature (SFT) values using various methods.

NO.	HM		MM		SRM		Proposed Method	Reference SFT
	OLR	QR	OLR	QR	OLR	QR
Data 1	75.871	79.372	74.713	77.936	77.362	86.914	80.999	80
Data 2	119.66	123.26	117.45	121.15	125.63	131.91	124.231	120
Data 3	20.328	18.981	20.003	21.302	22.132	19.499	20.327	20.25
Data 4	187.60	223.5	220.64	230.65	207.39	260.87	236.084	240
Data 5	98.449	100.07	97.048	98.257	99.858	106.16	100.133	101.111
Data 6	94.261	98.535	94.971	N/A	94.334	101.23	105.296	105

. It is obvious that the SFT values predicted by the proposed method are much better than others.

4. Applications in the Case of a Small Number of Data Points

The fitting ability of the proposed method has been tested in the cases of a set having a small number of BHT data points. Three different analytical methods with two different regression models mentioned in Section 2 were applied for comparison, the results were shown in

Table 5. Deviation percentages of the SFT estimates using different methods.

Method		Data 1								Data 2
		Data Number								Data Number
		3		4		5		(n)		3		4		5		(n)
HM	OLS	−11.756		−10.378		−8.765		−5.161		−6.808		−5.583		−4.608		−0.283
	QR	−5.456		−3.944		−3.058		−0.785		0.067		0.658		1.217		2.717
MM	OLS	−10.186		−9.576		−8.091		−6.609		−7.708		−6.708		−5.908		−2.125
	QR	−6.975		−5.304		−3.828		−2.580		−2.933		−2.100		−1.508		0.958
SRM	OLS	0.734		−7.094		−5.307		−3.298		−3.142		−1.517		−0.250		4.692
	QR	8.970		5.601		7.781		8.643		10.367		11.017		11.592		9.925
Proposed method		−2.367		0.028		0.470		1.249		3.929		3.971		4.289		3.527
Method		Data 3					Data 4						Data 5
		Data Number					Data Number						Data Number
		3	4		5	(n)	3		4	5	all		3	4	5		(n)
HM	OLS	1.600	0.835		0.212	0.385	−36.800		−31.438	−26.250	−21.833		−5.350	−4.717	−4.125		−2.632
	QR	1.215	−16.928		−10.731	−6.267	−14.917		−12.404	−9.500	−6.875		−4.779	−2.457	−1.191		−1.029
MM	OLS	−0.435	−3.289		−2.528	−1.220	−14.583		−11.029	−8.067	−8.067		−5.123	−4.740	−4.602		−4.017
	QR	0.652	2.588		2.074	5.195	−8.842		−6.358	−3.896	−3.896		−2.332	−2.546	−3.131		−2.822
SRM	OLS	9.852	9.931		9.151	9.294	−28.617		−22.992	−10.754	−13.588		0.613	−2.937	−2.070		−1.238
	QR	6.963	−22.089		−11.501	−3.709	14.696		12.213	9.367	8.696		3.907	−0.386	3.106		4.995
Proposed method		−4.958	−4.844		−0.543	0.380	9.501		0.784	−2.198	−1.632		−3.667	−2.661	−0.781		−0.967

. The first three, four, or five data points of each dataset were applied to estimate SFT in different analytical methods including the proposed method. The results are discussed in this section.

Figure 5 and Figure 6 present the SFT estimation results using different methods when the first three, four, or five data points or all the n data points in each set (i.e., DN 3, DN4, DN 5, and DN n) are used. The red full line in Figure 5a–f (plotted based on the dataset 1 to 6) represents the true value of SFT for each dataset, whereas the black dashed lines above and below the red line represent the SFT values with the DEV% of 5% and −5% of the true SFT value, respectively. The circles in Figure 5 and the diamonds in Figure 6 represent the results of SFT estimations using different methods. It should be mentioned that the red diamonds in Figure 6 represent the SFT estimates of our proposed method. The results of deviation percentages between each estimate with the reference SFT value are listed in Table 4 for reference. The SFT estimates discussed here were obtained using seven methods when the first three, four, or five points, or all data elements in each set, were applied, in the ranges [SFT −5% SFT to SFT +5% SFT], which is called “the acceptable range”.

By comparing the SFT estimates with the reference SFT values, it can be observed in Figure 5 that the SFT values estimated by the Horner method (HM) when OLR model is applied were typically underestimated in all datasets except for Data 3. Additionally, the aspheric heat-flow method (SRM) with the OLR model was likely to overestimate the SFT, as mentioned in [17]. As for the regression models, results showed that SFT values estimated by using an ordinary linear regression model were less than those by the quadratic regression model in Figure 5. This observation was also reported in other literature [16,17].

For Data 1 and 2 set, the proposed method is the only method in which all the SFT estimates lie in the acceptable range. It can be concluded that both the HM and MM method with the QR model are two good estimating methods, except for the proposed method from Figure 5a,b. However, the SFT estimates obtained using the two methods for only the first three points of Data 1 set were both out of the acceptable range. Estimation results for the Data 3 set using the three methods were all in the acceptable range. The three methods in this case involved HM with the OLR model, MM with the OLR model, and the proposed method. For this dataset use of the SRM methods with the OLR and QR models could not yield good results for the estimated SFT. For Data 4, application of DN 3 in all methods resulted into the SFT being out of the acceptable range. However, it should be noted that the proposed method was the only method in which the estimates laid in the acceptable range (using DN 4, DN 5 and DN n). All but one SFT estimated by the MM method with QR model on the DN n were unacceptable. For Data 5, the predicting results of the proposed method and other four methods all lied in the acceptable range. For Data 7 without the reference SFT value, the predicting ability of each method can still be seen in Figure 5f. The estimates obtained from HM method with OLR model and SRM method with the OLR model were uncertain, with a range of 63.77 °C and 99.74 °C, respectively. However, the proposed method, as well as the MM method with the OLR model, can still be used to estimate SFT reliably. Except for the SRM method with the QR model, using our proposed method, the estimates using other methods were likely to be underestimated when only the first three, four, or five points of the dataset were used, while the estimates by the SRM method with the QR model were more likely to be exaggerated. However, estimates obtained from the proposed method were likely to distribute but close on both sides of the reference SFT value (red full lines in Figure 5), which is also an indication of the advantage of the proposed method.

Based on the results of these statistical tests and the analysis above, the proposed method can be considered to estimate SFT accurately and reliably in the cases with few data.

Both data quantity and method will affect the accuracy of prediction, and some methods under specific conditions may have better results than the newly-proposed method, so the Theil inequality coefficient (TIC), an overall evaluation, was used to assess all of the data. The value of the TIC is in [0,1], and could be used to compare the results of different datasets. For each method, the closer the TIC is to zero, the more accurate that method is. The results of the TIC are depicted in Figure 7. It is to be observed that the new proposed method is the only one in which the TIC values for all datasets are less than 3%, indicating that the newly-proposed method cannot only estimate SFT using a few data points accurately, but also reliably.

Note that the SFT may vary with depth and rock type of the formation because the corresponding thermal stress varies with depth and rock type of the reservoir. This, however, has not been investigated in the current study and a research plan has been embarked to further explore the possible relationship between SFT and depth, as well as rock type.

5. Conclusions

(1)

A numerical method was proposed to estimate SFT from the BHT data and shut-in time. The unknown coefficients of the model were derived from the least squares fit by the particle swarm optimization (PSO) algorithm.

(2)

The estimation accuracy and fitting ability of the proposed method was verified using eight BHT datasets, including synthetic, geothermal, and petroleum field data. The deviation percentages are less than ±4% and the regression coefficient R² are greater than 0.98.

(3)

A comparison among different methods was conducted. The new proposed method could estimate SFT accurately and reliably, even by using a small number of BHT data points (the TIC values for all datasets are less than 3%). This might be used as a practical tool to predict SFT in both geothermal and oil wells.