Uncertainty Assessment of Corrected Bottom-Hole Temperatures Based on Monte Carlo Techniques

Abstract

Most temperature predictions for deep geothermal applications rely on correcting bottom-hole temperatures (BHTs) to undisturbed or static formation temperatures (SFTs). The data used for BHT correction are usually of low quality due to a lack of information and poor documentation, and the uncertainty of the corrected SFT is therefore unknown. It is supposed that the error within the input data exceeds the error due to the uncertainty of the different correction schemes. To verify this, we combined a global sensitivity study with Sobol indices of six easy-to-use conventional correction schemes of the BHT data set of the Bavarian Molasse Basin with an uncertainty study and developed a workflow that aims at presenting a valid error range of the corrected SFTs depending on the quality of their input data. The results give an indication of which of the investigated correction methods should be used depending on the input data, as well as show that the unknown error in the input parameters exceeds the error of the individual BHT correction methods as such. The developed a priori uncertainty-based BHT correction helps to provide a real estimate of the subsurface temperatures needed for geothermal prospecting and probabilistic risk assessment.

1. Introduction

Knowledge of the static formation temperature (SFT) is essential for both geothermal and hydrocarbon projects for optimal borehole design (drilling and completion, e.g., [1]), for exploring heat deposits, as a basis for interpreting geophysical loggings, e.g., [2], and for the calculation and correction of geothermal and hydraulic parameters [1]. It is of special importance for geothermal projects in low-enthalpy regions, since the SFT directly correlates with the amount of energy produced by a geothermal well, which is crucial for reaching the minimum production temperature [3]. Small differences in expected formation temperature therefore have large impacts on the estimated efficiency and economics of geothermal systems and, consequently, on exploration risk and insurability. Statistically validated temperature forecasts are consequently of high interest for geothermal prospecting and to reduce the risk of successful exploration. To obtain an approximation of subsurface temperatures, isothermal maps and temperature distribution models are used. The quality of these temperature predictions is influenced by the data density, the heterogeneity of data, including disturbed and undisturbed temperatures, and the applied inter- and extrapolation [4]. Disturbed temperatures represent data that are measured immediately or shortly after the drilling process and are therefore thermally affected by the drilling process, while undisturbed temperatures represent the actual rock temperature. The direct measurement of undisturbed temperatures is costly and time consuming, as well as challenging due to the unknown and proportionally long shut-in time that must be waited for. This is possible with a temperature logging tool, which is usually operated via a wireline, or with fiber optic temperature measurements, which are either wireline or permanently installed. Therefore, the SFT is rarely known and in most cases must be estimated from poor quality temperature data.

Most common sources of downhole temperature data are bottom-hole temperatures (BHTs) that are measured during or after drilling as a byproduct of geophysical logging runs and are usually documented at the lowest point of a logging run across a well section. Usually, the drilling fluid cools the formation, but at shallow depths at around 1000 m, the rock temperature is low and BHT measurements may therefore exceed the SFT, e.g., [5] (see Figure 1).

In the early stages of geothermal projects, BHTs are often the only downhole temperature information available as documented, e.g., by [6] for their Denmark data set. BHTs are well known to be highly influenced by mud circulation (temperature of the drilling mud, duration and pumping rate), the time passed after circulation stops, the geometry of the borehole, the geothermal gradient, and the thermal properties of the borehole fluid and the surrounding formation, e.g., [1].

An estimation of the SFT from BHT data is possible by analytical or numerical extrapolation of these temperature buildup during the shut-in period of the respective well. Reviews of the most prominent existing conventional correction methods highlighting their differences are given by [7,8,9]. Widely used are the Horner plot method [10] adapted by [10,11], the correction according to [12], as well as the generalized Kutasov–Eppelbaum method [13]. These methods represent BHT-correction procedures based on a line source approach that determine the equilibrium temperature as the intersection point with the temperature axis. A mathematical issue with these models is that they inadequately represent the borehole as they ignore the spatial dilatation (borehole radius). Other methods that are solved graphically are the spherical-radial heat flow model of [14] or the radial heat source Brennand method [15], which is recommended by several authors due to its accuracy [9,16,17]. The numerical and analytical methods following [18] are based on a cylindrical source model considering the spatial dimensioning (radius ≠ 0) [8,19]. By estimating the initial temperature disturbance and thermal diffusivity, one can also calculate the SFT for only one available BHT at fixed values for the borehole radius and thermal diffusivity, e.g., [20,21]. The latter method was refined by [22] by correlating the thermal disturbance with the borehole radius.

More advanced and more complex methods were introduced by several authors in recent years, such as a two-component model that respects the thermal interactions between drilling mud and formation due to their different thermal properties [23] or an artificial neural network approach to empirically derive the SFT [24]. Other available approaches include Wong-Loya’s rational polynomial functions [25], the stochastic method of particle swarm optimization algorithm [26], the multicomponent approaches on the basis of geochemical data [27,28], or a numerical wellbore simulator and numerical inversion of temperature logs during injection [29,30].

However, conventional methods are still widely used, e.g., [31], because they are easier to apply and require less information about the geophysical properties of the drilling fluid and the formation, as well as less computational effort than the more complex methods. Thus, they are often applied at wells with poor documentation where the accuracy of the information may be doubted or where the necessary parameters for the correction methods are often not known at all. This is particularly true for old hydrocarbon wells. In general, BHT data found in older logging reports and logging headers are not reliable, as the values often suggest incorrect documentation and rounding of relevant parameters [7]. A validation of different conventional BHT correction schemes versus field data was executed by, e.g., [8,32], who estimated errors up to about 8 to 10 °C, or [33] who studied different regression schemes for the solution of graphical BHT correction. The author of [32] concluded that, due to high standard deviations, corrected BHTs reflect that the undisturbed formation temperature is not better than ±16 °C. These authors focused on the different mathematical correction methods but did not include the likely interactions of varying input parameters in their studies. We suspect that the inaccuracy of the input data has a much greater impact on the correction than the methodological weakness of the generalized conventional methods. This was also highlighted by, e.g., [7], who stated that the accuracy of the corrected BHTs does not necessarily lie in the method used for the correction, but rather in the quality of the available data. The authors of [34] also found that the results of different correction methods differ a lot and that the errors are generally of importance, with up to 10 to 20 °C. They see the different mud circulation practices as a main trigger. For a holistic consideration of the true error of the corrected BHTs, the errors in all input parameters are relevant.

For prospecting and risk assessments prior to the start of new geothermal projects and for a general assessment of a region’s geothermal potential, tools that present reservoir temperatures as a reliable business case are needed. To date, a quantitative consideration of the error of BHT-corrected values has been lacking but is crucial for a risk assessment of the temperature forecast. Therefore, in this study, we investigate how BHT-corrected values can be provided as a distribution dependent on the quality of the input data and their errors, rather than as a distinct corrected temperature value that ignores error propagation, as is commonly done, e.g., [6,9,22,26]. To do so, we study the sensitivity of typically used correction schemes and their error considering different parameter assumptions. Based on the results of the sensitivity study, a new workflow that accounts for the likely errors in the input parameters and provides a probabilistic alternative to existing temperature estimates to assess business, as well as the worst- and best-case scenarios, is described.

2. Materials and Methods

We applied a variance-based global sensitivity analysis (GSA) using Sobol indices [35] and Saltelli sampling [36] on the temperature data set of the Bavarian Molasse Basin to show possible sources of error of BHT correction and provide recommendations on the favored methods to be applied. The focus was on six conventional BHT correction schemes that are commonly used and can be applied on a large data set with a partly poor or unknown quality, and thus large uncertainty of the available data. Therefore, the uncertainties in all input parameters were studied and the outcome’s uncertainty based on the given input parameter set and its error-proneness was evaluated. Using undisturbed temperature logs, also derived from fiber-optic distributed temperature sensing measurements (DTS) from the data set, the BHT correction schemes were compared with respect to their uncertainty. For this purpose, the possible deviations of the corrected SFT from the actual temperature conditions were calculated from the temperature logs and the DTS measurements, respectively. By combining the results of the sensitivity study and the uncertainty analysis, we developed a workflow that yields a corrected SFT with statistically valid uncertainty depending on the quality of the individual BHT data sets. Since conventional correction methods are only based on conductive heat transfer, we also investigated the performance of these BHT corrections in a well at known inflow zones.

2.1. Study Area and Data Set

The main exploration target in the conductive-dominated hydrogeothermal study area is the fractured and karstified carbonate rocks of the Upper Jurassic Malm aquifer. The production temperatures range between under 50 °C in the north, where the aquifer is at shallow depths, and over 160 °C in the south, where the Upper Jurassic is in depths of 4000 to 6000 m [37,38]. We used this temperature data set of drillings in the Bavarian Molasse Basin in Southern Germany, which includes both hydrocarbon and geothermal wells. In the 1950s and 1960s, multiple gas and oil wells were being drilled there, exposing layers of the Cretaceous and deeper in layers of the Upper Jurassic. In the last 15 years, a lot of geothermal projects were developed, making the Bavarian Molasse Basin one of Europe’s most important low enthalpy geothermal hotspots today [39,40]. Extensive studies were conducted during the geothermal buildup to characterize the reservoir of the Molasse Basin in detail as well as regionally, hydrochemically [41], petrophysically [42], and geothermally [43,44].

In the study area, a total of 65 geothermal wells and 870 hydrocarbon wells exist. Of those, BHT data with associated reported shut-in times from logging reports and headers are available from 346 wells (292 oil and gas wells and 54 geothermal wells). Fourteen wireline temperature logs (TLogs) of both geothermal and hydrocarbon wells were available, of which eleven were measured after a shut-in period of at least two months. In addition, a time series of continuously measured high-quality temperature profiles was available at the well TH4 of the Schäftlarnstraße (SLS) site in Munich. This well is equipped with a permanent fiber-optic monitoring system, installed inside of the borehole over its entire length along a steel sucker rod [45]. There, distributed temperature sensing (DTS) was used to collect temperature profiles over a 16-month period during the shut-in and to monitor the equilibrating borehole temperature. The DTS temperature profiles were acquired every 10 min at a spatial resolution of 1 m. With the exception of a dominant hydraulically active zone at the top of the reservoir, thermal dynamics (warm back) are no longer evident in the profiles after 16 months of shut-in [45]. Therefore, we suspect that this profile reflects the geothermal gradient with great accuracy. Further details on the background of downhole fiber-optic temperature measurements, the installation of the monitoring system in SLS TH4, and the technical specification can be found in [45]. A summary of all available continuous temperature logs is given in Table S1. Figure 2 depicts all available temperature logs (wireline and DTS) and the BHT data of the Bavarian Molasse Basin.

As can be seen in Figure 2, the disturbed BHT readings generally tend to overestimate the undisturbed temperature profiles at shallower depths and underestimate them at greater depths. For 165 cases, we had to assume copying errors of the BHT (i.e., the same BHT values were written down for different shut-in times), and in some other cases we noticed errors copying the shut-in time. These data were not used for further analysis. In total, there were 730 cases with one or a series of BHT measurements with reported shut-in times that could be corrected to SFT. For 155 wells, only a single BHT was available at a depth layer of the respective well. This means that more than one BHT may be measured in a well, e.g., from different logging runs at different depths/borehole sections of the respective well, but not consecutively in a series of logging runs in a borehole section. In the following, such data are then referred to as “one BHT.” Data for one BHT were generated if only one logging run was performed per borehole section or if a series of logging runs was performed but only one valid BHT was reported. If multiple consecutive logging runs were performed, but the reported BHT data did not increase over time, these BHT data were classified as invalid.

Table 1. The BHT data set structured in steps of 500 m TVD regarding the shut-in time.

Depth (mTVD)	Number of Wells	Valid Bottom-Hole Temperature at One Depth						Average Shut-in Time (s)	Temperature Range (°C)
		1 BHT	2 BHTs	3 BHTs	4 BHTs	5 BHTs	6 BHTs
≤500	39	36	7	0	0	0	0	30,600	21.0–61.3
500–1000	98	78	25	7	1	0	1	29,880	22.0–67.0
1000–1500	104	63	32	12	9	2	0	33,120	30.0–70.0
1500–2000	103	82	26	10	5	1	0	43,200	27.8–82.0
2000–2500	95	79	22	9	3	3	1	48,960	39.0–92.0
2500–3000	66	47	18	5	5	3	0	60,480	47.0–116.7
3000–3500	37	25	10	6	2	3	1	61,200	60.0–124.0
3500–4000	27	17	11	3	3	0	1	76,680	61.0–137.5
>4000	35	30	14	7	3	1	0	57,960	73.0–167.7

shows the evaluable BHT data in increments of 500 m TVD.

Even for wells drilled in recent years, the vast majority have only one BHT that can be corrected. A reason for this lies in the fact that often only a few expensive logging runs are performed and the BHT and corresponding shut-in time are reported only as a byproduct to the actual logging parameters and are thus often imprecise or missing.

In nine cases, a complementary data set of undisturbed temperature data was available from TLog or fiber optic DTS data (well SLS TH4) and BHT data with reported shut-in times. Unfortunately, only in two of these cases was a set of measurements from at least two BHTs at the same depth with different shut-in times available. The unperturbed temperature information was used as the target value in the uncertainty study in these nine cases.

2.2. Applied Correction Schemes

We applied the six conventional and commonly used correction techniques of the Horner Method (HM), the Lachenbruch and Brewer Method (LBM), the Brennand Method (BM), the linearization method (LM), forward modeling (FM), and the 1BHT Method (1BHTM). HM and LBM are based on a linear heat source, disregarding the radius of the borehole. These methods, as well as BM, based on a radial heat source, are solved graphically by interpolation of a time equivalent versus the measured BHT values. The LM, FM, and 1BHTM schemes, in contrast, are based on Leblanc’s cylindrical heat source model following Equation (1) [18]:

$$BHT \left(0,t\right)=SFT+\Delta T\left(e^{\frac{-a^2}{4\kappa t_s}}-1\right), \Delta T=SFT-T_m$$

(1)

where $BHT \left[°\mathrm{C}\right]$ represents the thermally disturbed in situ temperature, $SFT \left[°\mathrm{C}\right]$ the undisturbed rock temperature measured at equilibrium conditions, $\Delta T \left[K\right]$ the initial temperature disturbance, $a \left[m\right]$ the borehole radius, $\kappa \left[\frac{m^2}{s}\right]$ the bulk thermal diffusivity of the system drilling mud and formation, $t_s \left[s\right]$ the shut-in time (time since drilling fluid circulation stopped), and $T_m \left[°\mathrm{C}\right]$ the mud temperature (temperature of the drilling fluid during circulation). The analytical solution of Leblanc’s correction technique of a cylindrical explosion heat source was designed for at least three BHT data per depth. The method can be applied if the following stability criterion is met:

$$a^2<4\kappa t_s$$

(2)

The six applied methods are shown in Figure 3 subdivided by their mathematical approach.

We can also distinguish the six methods by how they are solved, i.e., analytically or graphically.

Table 2. Input requirements (marked with x) of six conventional BHT correction methods.

Graphically Solved
Name	Number of BHTs per Depth Layer	Shut-in Time	Mud Temperature	Borehole Radius	Circulation Time	Thermal Diffusivity
LBM	≥2	x
HM	≥2	x			x
BM	≥2	x			x
Analytically Solved
LM	≥3	x		x
FM	≥2	x		x		x
1BHTM	1	x	x	x		x

gives an overview of the six methods with their respective required input parameters.

2.2.1. Lachenbruch and Brewer Method (LBM)

The correction according to Lachenbruch and Brewer [12] determines the SFT as the intersection point of the BHTs with the temperature axis. Therefore, we used linear regression considering the reciprocal value of the shut-in time as $\ln \frac{1}{t_s}$ on the x-axis. For this method, at least two BHT values with associated shut-in times are required.

2.2.2. Horner Method (HM)

In contrast to LBM, the correction method according to Horner, originally developed for pressure build-up [10] and adapted for temperature build-up [11,46], additionally considers the circulation time $t_c \left[s\right]$ in the form of $\ln \frac{t_c+t_s}{t_s}$ on the axis of abscissa. In [48] it was shown that Horner’s method gives poor results for large borehole radii and short shut-in times, but reliable ones when the following criterion is met:

$$a^2\ll 4\kappa t_s$$

(3)

2.2.3. Brennand Method (BM)

The Brennand method is based on a radial heat source and is written as [15]:

$$\begin{gathered} BHT\left(a,t\right)=SFT-m\log \left(t_B\right), t_B=\frac{1}{t_s+p{t}_c}, \\ m=\beta \left(SFT-T_m\right)n, n=\frac{c_{p,r}{\rho }_r{a}^2}{\lambda } \end{gathered}$$

(4)

where $t_B \left[s\right]$ is the Brennand time, $p=0.785$ and $\beta =6.28$ are constants derived from field data of a Philippines data set, and $c_{p,r} \left[\frac{J}{kgK}\right]$, ${\rho }_r \left[\frac{kg}{m^3}\right]$, and $\lambda \left[\frac{W}{mK}\right]$ are the specific heat capacity, density, and thermal conductivity of the formation rock. The solution is derived when at least two BHTs recorded at different shut-in times are plotted versus their respective Brennand times. The intersection at zero Brennand time of a linear fit through all points represents the undisturbed formation temperature $\mathrm{SFT}$.

2.2.4. Linearization Method (LM)

LeBlanc’s formula (Equation (1)) can be solved inversely by numerical optimization of the thermal diffusivity [47,49]. Transformation and logarithmic calculus of Equation (1) generate the linear equation according to the scheme C = a − b ∗ t*:

$$\underset{C}{\underbrace{\ln \left(BHT-T_m\right)}}=\underset{a}{\underbrace{\ln \left(\Delta T\right)}}-\underset{b}{\underbrace{\left(\frac{a^2}{4}\right)\ast \frac{1}{\kappa }}}\ast \underset{t*}{\underbrace{\frac{1}{t_s}}}$$

(5)

where SFT = Tm + ∆T. Equation (5) can be transformed into an inverse linear optimization problem in case at least three independent BHT measurements are available. In that case, initial guesses need to be made for the mud temperature $T_m$ and the bulk thermal diffusivity of the drilling mud and the nearby rocks $\kappa$, while the borehole radius a, the shut-in time $t_s$ and the measured BHT values are given. As a first approximation, $T_m$ can be set at 30 °C or at 50% of the maximum observed BHT value.

During the iterative linear optimization procedure applying standard LSQ methods, the correction vector adapts $\Delta T$ and $\kappa$ until a convergent solution is found or a maximum number of iterations has elapsed. The resulting value for $\kappa$ can be used for quality control measures to evaluate if thermal conduction is the major process of the thermal balancing observed inside the borehole regarding the observed BHTs. In addition, this method delivers the total fitting error (predicted versus observed BHT values) for all iteration steps.

2.2.5. Forward Modeling (FM)

In contrast to the linearization method, forward modelling is an analytical solution that disregards mud temperature. The method can be used for at least two BHTs per depth. The following Equation (6) represents the calculation of undisturbed rock temperatures by analyzing the temperature increase with raising shut-in times [20]:

$$SFT=BHT\left(t_{s,1}\right)+\left(BHT\left(t_{s,2}\right)-BHT\left(t_{s,1}\right)\right)\frac{1-f\left(t_{s,1}\right)}{f\left(t_{s,2}\right)-f\left(t_{s,1}\right)}, f\left(t_s\right)=e^{\left(-\frac{a^2}{4\kappa t_s}\right)}$$

(6)

2.2.6. 1BHT Correction Scheme (1BHTM)

All the methods described above have in common that they are applicable only when at least two BHTs have been measured at the same depth and at different time intervals. In most boreholes (see Table 2), only one BHT is available at a depth layer. Leblanc’s Equation (1) can be simplified to also account for the correction of a single BHT. The transformed equation writes as follows [18]:

$$SFT=\frac{BHT+T_m\left(f\left(t_s\right)-1\right)}{f\left(t_s\right)}, f\left(t_s\right)=e^{\left(-\frac{a^2}{4\kappa t_s}\right)}$$

(7)

2.2.7. Constraints

There are two restrictions to the correction schemes. First, a basic requirement for all methods is that each consecutive temperature measurement and the corresponding time since circulation stopped is higher than the previous one. Second, as some of the applied methods are valid only if the stability criteria $a^2<4\kappa t_s$ (Equation (1)) or $a^2\ll 4\kappa t_s$ (Equation (3)) are fulfilled, we applied the respective criterion to all methods to keep the results of the sensitivity study comparable. Fulfillment of the stability criterion is dependent on borehole radius $a$, shut-in time $t_s$, and estimated thermal diffusivity $\kappa$. We therefore applied the sensitivity study at typical borehole radii, as they are typical for the Bavarian Molasse Basin, and calculated the limits of minimum and maximum shut-in time and thermal diffusivity with respect to the criterion. The typical diameters of the well sections in the data set were 23 inches (0.58 m), 17.5 inches (0.44 m), 12.25 inches (0.31 m), 8.5 inches (0.22 m), 6.25 inches (0.16 m), and 6 inches (0.15 m). For each section, we considered a possible widening of 1 inch (~0.03 m) due to outcrops in the rock.

2.3. Global Sensitivity Analysis (GSA)

The quality of input parameters to the BHT correction methods within the data set varied. For example, even-numbered shut-in times without a decimal digit frequently occurred in logging headers, which indicated to us that these values were often rounded or not measured precisely. In other cases, the values were more accurate to one decimal place, implying that they were more reliable. We applied a global sensitivity analysis (GSA) using Sobol method with the open source Python library SALib [50] and Saltelli sampling [36,51] to the BHT correction methods (HM, LBM, BM, LM, FM, and 1BHTM) to investigate the qualitative influence of parameter assumptions or empirical approaches. Sobol method is a robust and high-performance [52] variance-based GSA for which all input parameters are varied over the whole parameter space [35]. The analysis produces first-order indices that determine the impact (percentage) of the variance of an input parameter on the model output, second-order indices that measure the interaction between two parameter variances, and total-order indices that determine the overall effect, including interactions that a parameter variance has on the entire model output. The suitability and strength of GSA for geoscience applications with the Python library SALib has been demonstrated in several studies, e.g., [53,54]. To formulate the problem for Sobol analysis, we defined uncertainty ranges for each input parameter set, then chose a realistic statistical distribution for each set and sampled it according to Saltelli’s extension of the Sobol sequence [36,51]. SALiB makes it possible to specify the parameter sets as four basic distributions, which are rectangular, if all parameters within the set are equally likely, or triangular, normal, or lognormal, if the parameter is non-uniformly distributed.

2.3.1. Parameterization of Borehole Radius a

For the borehole radius, a triangular distribution was assumed, with the minimum and mean value as the borehole radius resulting from the respective well section diameter and the maximum value being an additional 0.03 m, which is possible due to outcrops but rare in our study area.

2.3.2. Parameterization of Thermal Diffusivity κ

The applied correction methods FM and 1BHTM require the thermal diffusivity κ as a bulk value for the whole borehole system (drilling mud/formation). κ is dependent on the specific geological setting and cannot be measured downhole. Thus, it must be estimated. To choose a realistic min/max range, we researched literature values for the bulk thermal diffusivity of carbonate and sedimentary rocks as found in our study area of the region of the Bavarian Molasse Basin. Within those, we found a minimum value of 1.5 × 10⁻⁷ m²/s assumed from numerical tests and statistical data in [22] and a maximum value of 6.8 × 10⁻⁷ m²/s in [48]. Other values we researched from [5,18,55,56,57] lie in between these extrema.

We varied κ between 1.5 × 10⁻⁷ m²/s and 6.8 × 10⁻⁷ m²/s in agreement with the researched literature values and distributed them uniformly by a rectangular type.

2.3.3. Parameterization of Shut-in Time t_s

The shut-in time $t_s$ is a required parameter for any correction scheme. From the selected range for κ, we calculated the minimum shut-in times for which the stability criteria Equation (2) is still met. This is true if the squared radius $a^2$ in the left term of Equation (2) is smaller than four times the product of $t_s$ and κ in the right term of Equation (2). The calculated minimum shut-in times are shown in

Table 3. Calculated minimum shut-in times for which the stability criterion is still met when κ is varied between 1.5 × 10⁻⁷ and 6.8 × 10⁻⁷ m²/s.

Φ Bit Size [inch/m]	a2 [m2] (with Additional 1 Inch Due to Outcrops)	${\mathit{t}}_{\mathit{m}\mathit{i}\mathit{n}}$ [s] for ${\mathit{\kappa }}_{\mathit{m}\mathit{i}\mathit{n}}=1.5\times {10}^{-7} {\mathbf{m}}^2/\mathbf{s}$	${\mathit{t}}_{\mathit{m}\mathit{i}\mathit{n}}$ [s] for ${\mathit{\kappa }}_{\mathit{m}\mathit{a}\mathit{x}}=6.8\times {10}^{-7} {\mathbf{m}}^2/\mathbf{s}$
23.00/0.58	0.1008	168,010	37,061
17.50/0.44	0.0613	102,218	22,548
12.25/0.31	0.0328	54,587	12,041
8.50/0.22	0.0178	29,637	6538
6.25/0.16	0.0110	18,296	4036
6.00/0.15	0.0103	17,204	3795

for the extrema of the κ range.

Figure 4 shows that the data set is represented well at the smaller borehole sections (8.5 inches, 6.25 inches, and 6 inches) for which the majority of reported shut-in times exceeds the calculated minimum shut-in time at minimum κ. For 12.25 inches 25%, for 17.5 inches only 3%, and for 23 inches no BHT data are represented.

The maximum shut-in time of all BHTs recorded in the data set was 170,000 s, and the minimum and maximum hold-up times between a previous and a subsequent recording were 1800 s and 125,000 s, respectively. For GSA, we therefore varied the first recorded shut-in time between each minimum value of Table 3 (depending on the borehole radius used) and 170,000 s. For every subsequent measurement, the respective shut-in time must be larger than the previous. To represent this in the GSA, we filtered the data set for cases showing more than one BHT per depth and calculated the deltas and implemented the second and third values as a positive delta (delta_t_i) to the respective preceding value. The range for delta_t1 and delta_t2 were subsequently set to the minimum waiting time from the previous measure of 1800 s and maximum waiting time of 125,000 s and distributed uniformly.

2.3.4. Parameterization of Measured In Situ Temperature BHT

The BHT value is a required input for any correction scheme. Although the BHT was correctly and accurately reported in the field after the geophysical logging run was performed, it is still subject to a measurement error due to the uncertainty of the tool. The uncertainties of common temperature gauges used at logging tools are in the range of ±2 °C [58] or ±3% (kind note of the service company Weatherford). For the GSA, we chose each initial BHT as a representative value for each well section from the data set (6 inches: 90 °C, 6.25 inches: 70 °C, 8.5 inches: 66 °C, 12.25 inches: 55 °C, 17.5 inches: 49 °C), and set the ranges for the input parameters in a conservative manner with an assumed error of ±3%. We distributed the BHT uniformly as a systematic measurement error is equally likely within the error range. As in a row of measurements subsequent BHTs are expected to be higher than those previously measured, we had to reflect this in the GSA. Consequently, we implemented subsequent BHTs as a positive delta (delta_iBHT) to the respective preceding BHT with a range of 0.5 K to 36 K, which was derived from the data set after filtering for cases showing more than one BHT per depth.

2.3.5. Parameterization of Circulation Time t_c

The circulation time $t_c$ is a required parameter for HM and BM. For the given data set of the Bavarian Molasse Basin (see Figure 2 and Table 1), the circulation time has been reported in very rare cases. From the few drilling reports on hand, it appears that the average minimum circulation time is about two hours and increases with the depth of a well. A relation of depth and circulation time seems to be recognizable, since we know that a higher volume requires a longer cleaning period and in general a longer subterranean intervention, which increases the thermal disturbance. In previous studies, an estimation in the form of Equation (8) was applied, assuming an increase in $t_c$ by two hours per 1000 m depth [59]:

$$t_c \left[s\right]=depth \left[m MD\right]\ast 7200 \left[s\right]$$

(8)

However, in order to represent a broad spectrum of possible circulation times in the GSA, we initially chose a wide range with 3600 s and 144,000 s as extrema and 7200 s as the most frequent value known from drilling reports as the peak of a triangular distribution.

2.3.6. Parameterization of Mud Temperature T_m

The mud temperature $T_m$ is one of the input parameters for the LM and 1BHT-correction scheme. In rare cases, mud reports that provided information about the inlet and outlet temperatures of the drilling mud, the respective drilling depths, and the pumping rates, were available. As Figure 5 shows, a linear regression is in our case not suitable for the prognosis of the mud temperature, as only little data are represented by the regression.

The mean of all mud outflow temperature values at hand is 54 °C. As the heat loss or gain between drilling depth and the surface is unknown, an estimation of the mud temperature from the outflow temperature is prone to unknown errors. Therefore, we applied a broad range of 24 to 80 °C, covering the known drilling fluid outflow temperatures from the data set and implemented the distribution of the parameter in a triangular form with 54 °C as the peak value and 24 and 80 °C as the minimum and maximum values, respectively.

2.3.7. Summary of the Variation of the Input Data

The resulting variances of the inputs built the basis of the parameter space for the Sobol analysis.

Table 4. Applied ranges for all input parameters for the uncertainty analysis. The last column lists the respective input parameters with their individual designation.

Parameter	Applied Uncertainty	Variance Range as [;]	Labeling for Analysis
BHT value	Accuracy of measurement tool and rounding error	[BHT × 0.97; BHT × 1.03]	1st_BHT, delta_1stBHT, delta_2ndBHT
First shut-in time	Minimum: shut-in time for which stability criteria are still fulfilled, Maximum: from data set	[minimum time; 170,000s]	t1
Subsequent shut-in times	From data set	[1800 s; 125,000 s]	delta_t1, delta_t2
Radius	+1 inch possible due to possible outcrops	[radius; radius + 0.0254 m]	radius
Cir circulation time	Unknown, parametrization over whole parameter space of data set (drilling reports at hand)	[3600 s; 144,000 s]	tc
Mud temperature	Unknown, parametrization over whole parameter space of data set (mud reports at hand)	[24 °C; 80 °C]	Tm
Thermal diffusivity	Unknown, parametrization over a wide range from literature	[1.5 × 10−7 m2/s; 6.8 × 10−7 m2/s]	kappa

summarizes the origin of the uncertainty of the different input parameters and their designation in the model.

2.3.8. Sampling and Model Convergence

The convergence of the GSA model solution had to be proven after every Sobol model run. For this, we repeated the Sobol analysis to increase sample numbers until we found a stable solution (e.g., where the Sobol total order index did not change anymore as depicted in Figure 6). The required number of samples to achieve convergence varied between 10,000 and 90,000 for our models, depending on the number of input parameters, e.g., [50].

2.4. Uncertainty Study

To quantify the uncertainty for each method, wells were chosen for which both the target value (the SFT from available continuous TLog, respective DTS data from SLS TH4) and the quality of the input data set were well known. The available TLogs were measured after a minimum of 2 months of no circulation. With such long shut-in times, we assumed that temperature equilibration was at an advanced stage and the logs reflected the approximate undisturbed formation temperature in the well. We developed the Python Script BHT_Unct that contains the introduced correction schemes HM, BM, LBM, FM, LM, and 1BHTM and uses Saltelli’s extension of Sobol sequence sampling to create the parameter space. We proceeded with a parameterization of the input parameters similar to that for the GSA (Table 4), except in cases where additional or more detailed information was available. For example, if there were drilling reports available, we took the circulation time from there with an estimated uncertainty of 3600 s to respect rounding and imprecise reporting. If there were mud reports at hand, the mud temperature was estimated from the reported outflow temperatures. For shut-in time, we took the reported values from the respective BHT measurements and rated the accuracy of the reported values. In the absence of a decimal place, it was assumed that the value was imprecisely documented. Then, an uncertainty of ±7200 s was applied to make a conservative estimate and to account for possible rounding errors. Shut-in times reported with a decimal place seemed to be measured with more caution, and we respected this by specifying a higher quality with a lower uncertainty of ±3600 s. This was the case for well no. 3 and no. 7.

The criterion that consecutive BHT measurements increase in temperature and shut-in time had to be met. After the uncertainty ranges are created, the parameter spaces may overlap, which may result in the parameter space sampled by Saltelli not satisfying this constraint. In this instance, the overlap region was therefore equally trimmed from both sides.

After we ran the different models, we studied the distribution of the model result space to find the uncertainty as the deviation from the known SFT. The percentiles for which 10%, 50%, or 90% of the data lie within (p10, p50, and p90 limits) were chosen to describe the distribution of the values in a probabilistic way. In this context, the p50 limit describes the median of all calculated values in the distribution and can be used as the expected value. The p10 and p90 values can be used as a worst-case and best-case prediction, respectively.

We examined seven wells (well no. 1 to well no. 7) with the known SFT from TLogs and DTS for 1BHTM, for which only one BHT value was measured, and the stability criterion Equation (2) was met. An overview of the wells, and the available BHT and SFT data, are given in Table S1. Figure 7 shows the data set of the well SLS TH4 with fiber-optic temperature data (16 months of shut-in well) and one BHT (105.4 °C) with a reported shut-in time of 86,400 s available at the bottom of the reservoir section.

From these data, it is concluded that the SFT at the respective depth is 109.4 °C. To provide a more robust analysis, BHT data from a nearby well at the Schäftlarnstraße site was also included. To investigate the uncertainty of methods that require more than one BHT per depth, the available data limit the analysis to two wells (no. 8 and no. 9, see Table S2) for which a long-term temperature log over the shut-in period is available and for which the reported shut-in time is large enough to meet the stability criterion. The temperature log of well no. 8 was measured two years after the last circulation; for well no. 9, a temperature log of a nearby well at 2300 m in distance was used to be measured after 2 months of shut-in. For well no. 8, the two BHTs were measured at 2240 m depth, where the SFT from TLog is assumed at 66.50 °C; for well no. 9, two BHTs were available at 1492 m MD depth, where an SFT of 81.25 °C was derived from the nearby TLog.

The application of the LM is restricted to these examples as it is designed for at least three BHTs at one depth. Therefore, a series of four BHTs at 2355 m MD in well no. 8 was used to analyze the uncertainty of LM. The TLog available in well no. 8 is incomplete and does not cover the depth at which the four BHTs were measured. Without known SFTs, the uncertainty of the LM method cannot be quantified, but it can be estimated in a qualitative way by comparison with other correction methods applied to the four BHT series (e.g., FM and BM that have been studied before at well no. 9 and the two available BHTs in well no. 8 at 2240 m depth and known SFT).

2.5. BHT Correction at Flow Zones

We studied correction methods that are based on conductive heat transport process only. It is well known that, for example, Horner’s method fails in the presence of strong convective processes, since the formation then takes longer to reach complete temperature equilibrium [6,14,60]. To conclusively assess the quality of a BHT correction using conventional methods, one should consider convective zones in the wellbore to ensure that the corrected BHT is not within one. To do so, the fiber-optically monitored well SLS TH4, which has a highly active hydraulic zone between 2820 m MD and 2900 m MD [45], was studied. We took DTS profiles that were measured after a 24 h lasting cold-water injection and applied Horner’s and Brennand’s method on temperatures at different depth intervals from different times after injection (10,800 s to about 2,340,000 s shut-in time). As Figure 8 shows, the corrected values above 2800 m MD and below 2950 m MD are close to the undisturbed DTS profile but have a high deviation (up to 4.2 K) from the undisturbed DTS profile in the known convectively dominated zone.

3. Results

We applied the GSA at the borehole sections 6 inches, 6.25 inches, 8.5 inches, and 12.25 inches, as the underlying parametrization is not representative for our data set at the 17.5-inch section and the 23-inch section (see Figure 4). We completed the GSA for all six methods for up to three available BHTs per depth. The total order indices are shown in Figure 9, classified by the graphical methods (HM, BM, and LBM) and the analytical methods (LM, 1BHTM, and FM).

For all graphical methods, the circulation time has low sensitivity (total order index <0.1). In general, the sensitivity of shut-in time is higher (the total order index of second and third recorded value is between 0.25 and 0.5 for HM, BM, and LBM) than the BHT value (the total order index is maximum 0.3 for the second BHT measured for HM). For the 6-inch, 6.25-inch, and 8.5-inch radii, the results are very similar, whereas for 12.25 inches, the first measured shut-in time appears to be less sensitive than for smaller borehole radii. In general, the chosen range/uncertainty of the shut-in times has a larger impact than that of the BHT values (especially the second and third measured shut-in times, which have a total order index for HM and BM between 0.25 and 0.5, and the first measured shut-in time for LBM, which has a total order index of about 0.6).

For all analytical methods, the borehole radius has low sensitivity (total order index <0.1). In general, the FM tends to be more sensitive to the shut-in time (especially the third measured one with a total order index up to 0.7) at the chosen range/uncertainty of the input, whereas the LM tends to be more sensitive, especially to the second measured BHT value (total order index up to 0.62). The outcomes of both methods are stable for each section where the GSA was applied. In contrast, we see a varying sensitivity of 1BHTM in the different sections. The sensitivity of the mud temperature is lower at the 6-inch section (total order ~0.25) and higher at larger radii (e.g., at 12.25 inches total order ~0.4). The reverse is true for the measured shut-in time (total order > 0.5 at 6 inches and ~0.25 at 12.25 inches). The sensitivity of the thermal diffusivity is between the total order 0.25 and the total order 0.38 for all sections. Compared with those three input parameters, the BHT as an input value has less influence on the result (total order < 0.25 at all sections).

3.1. Uncertainty of 1BHTM

The input parameters of 1BHTM were the measured BHT, the shut-in time, the borehole radius, the thermal diffusivity, and the mud temperature. Given the uncertainties that must be assumed for the individual inputs, we can describe the total uncertainty of each solution as a deviation from the static formation temperature. The resulting ranges of the Sobol uncertainty analysis are shown as density and box plots in Figure 10 for the seven observed wells at different variances for the thermal diffusivity κ.

For the first plot from the left in Figure 10, κ was varied as in the GSA (1.5 × 10⁻⁷ m²/s to 6.8 × 10⁻⁷ m²/s); in the second and third plots, κ was set to different base values and varied by 50%. In all the observed cases, the SFT lies in between the p10 and p90 limits of the results. The maximum deviation of the result at base values from the respective SFT for κ = 3.0 × 10⁻⁷ m²/s is 8.3 K for well no. 6 and the minimum deviation is less than 0.5 K for well no. 4. The largest range between p10 and p90 limit was modelled at 31 K for well no. 5.

3.2. Uncertainty of Corrections for at Least Two BHTs at One Depth

To investigate the uncertainty of the HM, BM, LBM, and FM correction schemes, we performed model runs for two geothermal wells (no. 8 and no. 9), for each of which two BHTs and an undisturbed TLogs covering the depth of the measured BHTs were available. Figure 11 and Figure 12 show box diagrams and density plots of the result spaces for each model with its p10, p50, and p90 limits (first, second, and third black dashed lines), its modal value (gray line), the respective SFT (red dashed line), and the calculated value when the particular model is run without varying the input parameters (blue dashed line).

From the graphical methods shown in (a), (c), and (d) of Figure 11 and Figure 12, it can be observed that BM outperforms HM and LBM as the p50 value has a smaller deviation of 1.4 K (well no. 8), respective to 2.05 K (well no. 9) from the SFT compared with 3.1 K and 3.5 K for LBM and HM. The range between p10 and p90, however, is larger at 10.3 K (well no. 8) and 13.6 K (well no. 9) for BM compared to 7.5 K (well no. 8) and 11.8 K (well no. 9) for HM and LBM.

The ranges between p10 and p90 for FM are in the same order of magnitude as for BM, namely about 10 K at well no. 8 and about 14 K at well no. 9. The deviation of p50 from the SFT for FM is also in the same magnitude as for BM, namely about 1.5–2.0 K at well no. 8 and 2.1 K at well no. 9.

The linearization method LM was tested in well no. 8 on a BHT data set with four independent measurements, in comparison to the BM and FM correction schemes, as shown in Figure 13.

Since the four BHTs were taken at a depth where there was no TLog, and as this extended to only 150 m above the measured values, the accuracy of the LM method was evaluated in comparison to the FM and BM methods. Their uncertainty was previously investigated in the same borehole using a series of measurements with two BHTs (Figure 11 and Figure 12). LM appears to have the closest p10–p90 range of the three methods with 5.4 K in comparison to 5.5 K for BM and 6.7 K for FM. The calculated base value is of the same order of magnitude for all three methods; however, BM and FM seem to tend toward a higher deviation of the value at the p50 limit from the base value.

4. Discussion

By combining the findings of the GSA and of the uncertainty analyses, we can make recommendations for temperature prediction from BHT data. From the sensitivity study (Figure 9), it can be concluded that for the Bavarian Molasse Basin, the choice of the correction method to be applied should be made depending on the quality of the input data set of the individual wells. If the reported shut-in time is assumed to be reliable, both the forward modeling FM and Brennand method BM work well. If the reported shut-in time is likely to have a high uncertainty, it is recommended to use the linearization method LM if at least three valid BHT values are available. The 1BHTM correction scheme should be used only if there is not more than one consecutive BHT available. Figure 14 shows a proposed decision-making method regarding which method to apply depending on the number of available consecutively measured BHT values and on the uncertainty of the given input parameters.

The correction schemes used only work under conductive conditions (see Figure 8). Thus, we suggest that corrected BHT values measured in the reservoir at depths where flow zones are suspected or interpreted should not be corrected using conductive BHT correction schemes.

After a method is selected according to Figure 14, the BHT values can be corrected with the Sobol method for the given uncertainties of the input parameters. This procedure was implemented into the Python tool BHT_Unct, which is based on the Python library SALiB and can be obtained via GitHub (https://github.com/Flix-S/BHT_Unct, accessed on 15 July 2022) under an open-source GPL-3.0 license. The advantage of this approach is that the corrected values can be represented as a distribution function (as in Figure 10, Figure 11, Figure 12 and Figure 13), so that the statistical values of the density plot can be used to describe a business case (p50 of density plot), a worst-case prediction (p10), or a best-case prediction (p90).

4.1. Sensitivities of Parameters and Correction Schemes

The results of the GSA (Figure 9) show that one should evaluate the data set for the accuracy of the input parameters before choosing a method to correct BHT to SFT. The different input parameters contribute to the model results as follows.

In Figure 9, we showed that the shut-in times are highly sensitive in all graphical methods (HM, BM, and LBM). For HM, this is in accordance with the findings of, e.g., [23] or [48], which state that the method becomes more precise the longer the waiting time between shut-in and the BHT measurement. In some studies, the shut-in time was estimated using linear regression if no time was reported for the respective measurement, e.g., [22] (based on [61]). This seems reasonable, since longer shut-in times can be assumed for greater depths and higher temperatures, because the duration for extracting the drilling tool and then retracting the measurement tool increases correspondingly. However, Figure 15 shows that a linear regression based on depth, i.e., BHT measurement, does not represent the data set of the Bavarian Molasse Basin.

Given this fact and the high sensitivity, the use of the graphical correction schemes is not recommended if the reported shut-in time has a high uncertainty or if it is unknown.

The circulation time, which is required as input for HM and BM, has a marginal relevance with a total order index <0.1 (see Figure 9). Therefore, a rough estimation of this parameter following Equation (8) seems acceptable. Regarding the measured in situ temperature, we showed that the second measured BHT value is more sensitive to the graphical solutions than is the first measured one. This is explained by the fact that the second value significantly influences the slope of the regression line from which the SFT is estimated (see Figure 9).

For the analytical methods LM, 1BHTM, and FM, we can see a more diverse spread of sensitivities. The 1BHTM scheme is sensitive to $t_s$, $T_m,$ and $\kappa$. The borehole radius and BHT value are of minor sensitivity. As the sensitive parameters $T_m$ and $\kappa$ are usually unknown and hard to predict, the 1BHT method should generally be avoided in cases where other methods could be applied. $T_m$ must be estimated accurately, especially for corrections at larger radii. For FM, the shut-in times are important, especially the shut-in time of the latter measured BHTs. For LM, shut-in times have a minor role, but the measured BHT value, especially the second measured one, is sensitive to the model output. As an error of 3% can always be assumed in temperature measurement (see Section 2.3.3), it is recommended to use FM or BM and, if there are doubts about the accuracy of the reported shut-in times, LM.

4.2. Uncertainty

For the seven wells at which the 1BHT-method was applied, we found that an underestimation of the thermal diffusivity leads to an overestimation of the SFT (Figure 10). An excessively low assumed thermal diffusivity means a slower assumed spread of the temperature wave during drilling. This results in too high a correction. For the exemplary studied wells in the Bavarian Molasse Basin, 3.0 × 10⁻⁷ m²/s had the best fit for the thermal diffusivity. At this value, the mean of all Sobol results seems to be the estimation with the lowest uncertainty of the SFT.

The methods using two or more BHTs were tested at well no. 8 and well no. 9. In general, the maximum error of the graphical methods is high. These methods calculate the intersection of a regression line through the BHT values and the respective calculated representative time (Brennand time, Horner time, Lachenbruch time) at zero (BM, LBM) or one (HM). When two shut-in times or BHT values are varied up to their maximum assumed error (e.g., BHT1 − 3% error and BHT2 + 3% error), the resulting slope of the regression line changes on large scales.

For the two studied wells, we can support the statements that the Brennand method performs well [9,16,17] as it outperforms HM and LBM (see Figure 11 and Figure 12). Based on the results of only the two wells, it is not possible to deduce if BM or the analytical method FM should be preferred. At well no. 8, the deviation of p50 from the SFT is 2 K for FM and 1.4 K for BM, and at well no. 9, the deviations of FM and BM are equivalent, at about 2 °C.

A deeper BHT data set in well no. 8, for which no SFT could be obtained, was used to study the uncertainty of LM in comparison to FM and BM at the same data set. The distribution of results from FM and BM is similar, with a p50 value of 72.4 °C and 72.6 °C, respectively, and an 80% uncertainty range at 6.7 K and 5.5 K (Figure 13). LM, in comparison, has a similar uncertainty range of 5.4 K, but its p50 value is lower by 1.6 °C and 1.8 °C, respectively (Figure 13).

From the uncertainty analysis, we can conclude that very large errors are possible when an unfavorable combination of input parameter variances occurs. This is most evident in the tailing of the density plots, which are up to 60 K for 1BHTM (Figure 10, well no. 2, BHT no. 1) and about 30 K for BM at well no. 8 (Figure 11), for example. Such situations, however, should be recognized in the field if the available data is evaluated with caution, e.g., by filtering for temperature data or reported times that seem unrealistic for the completed logging runs or the known geological conditions at the site. Therefore, we expressed the likely error of the correction as the deviation of the most frequent value (modal value) or mean or p50 case from p10 and p90 cases, which should cover all realistic combinations.

4.3. Implications for Temperature Predictions in the Bavarian Molasse Basin

Only two TLogs were available at wells with more than one reported BHT per depth, and there was an incomplete data set of a well in the Molasse Basin with four BHTs at one depth and a known SFT from TLog or DTS. This emphasizes the importance of the 1BHTM correction scheme for this data set. The workflow used in this study proposes to use the p50 value of the distribution of corrected temperature as a business case for temperature predictions, p10 as a worst-case, and p90 as a best-case scenario. This procedure means that the corrected SFT at p50 can be specified at 80% with an uncertainty equal to the deviation at p10 and p90.

Thus, for the wells we studied with 1BHTM at best fit for $\kappa$ (3.0 × 10⁻⁷ m²/s, see Figure 10), 80% of the corrected values are within a range of 29.2 K maximum (well no. 2, BHT no. 1) and 7.6 K minimum (well no. 7). An example illustrates the relevance of these ranges of error. Assuming that the formation temperature is equal to the temperature of the produced water and that there are no heat losses along the production well, the SFT equals the production temperature $T_p$. With an exemplary pumping rate of $Q=0.1 {\mathrm{m}}^3/\mathrm{s}$, a fixed injection temperature $T_i$ assumed to be 50 °C, and an assumed heat capacity $c_f$ and density ${\rho }_f$ of the fluid of 4181 J/(kg*K) and 998 kg/m³, the significance of the temperature uncertainty for the heat output $P \left[W\right]$ can be calculated according to, e.g., [3]:

$$P=Q\ast {\rho }_f\ast c_f\ast \left(T_p-T_i\right)$$

(9)

If the best-case SFT prediction (p90) is assumed to be 100 °C, the output according to Equation (9) is 20.9 MW. For a maximum uncertainty range of 29.2 K, the worst-case scenario (SFT p10 at 70.8 °C) is then 8.7 MW. This means a clear reduction in thermal power by 58.4% if the worst case occurs, instead of the best case. For the minimum observed uncertainty range at well no. 7, the reduction in thermal power is calculated analogously to 15.2%. The ranges for HM and LBM are lower in wells no. 8 and no. 9, at 7.5 K and 11.8 K, respectively, representing a 15.0% and 23.6% reduction in thermal power when the same assumptions are made for Equation (9) as before. The 80% uncertainty ranges of BM and FM are in the same order of magnitude at about 10 K (20% output reduction) for well no. 8 and 13.6 K (27.2% output reduction), respectively. For the lower BHT data set of well no. 8, the uncertainty ranges of LM, BM, and FM are also in the same order of magnitude of about 5 to 7 K, implying a reduction in thermal output of about 10 to 14%. Such scales are of clear importance for economic and planning aspects.

However, comparing the methods is not possible by the 80% ranges alone. Figure 16 examined how the 80% ranges for the different methods change when the input parameters are assumed to have different uncertainty.

For the calculated cases on the right side of Figure 16, we assumed that the shut-in time was reliably documented and that an uncertainty of 900 s was applied to it (shut-in time ±900 s). The plots on the left show the calculated density plots for a higher uncertainty (shut-in time ±7200 s). The 80% range for 1BHTM (Figure 16 bottom) remains almost unchanged (15.4 K at low confidence of shut-in time and 14.5 K at high confidence), while the 80% range of FM clearly decreases by about 7 K. This shows the overall high uncertainty of 1BHTM, even if the quality of the input data set is satisfactory. On the other hand, it also shows the applicability of 1BHTM for BHT data sets of low quality, since the uncertainty of the other methods is then not significantly higher.

In contrast, if there is low confidence in the reported shut-in time, it is recommended to use the linearization method LM. As the circulation time is not sensitive to the studied models HM and BM, we propose to use the simple estimation of Equation (8) with the approach that the circulation time should last longer with rising depth. The thermal diffusivity and mud temperature are important input parameters when applying the 1BHTM. By comparing corrected BHTs with undisturbed TLogs and DTS in the Bavarian Molasse Basin (Figure 10), we found a κ of 3.0 × 10⁻⁷ m²/s as a best fit.

5. Conclusions

The results obtained confirm the findings of former studies according to which the accurate correction of BHT data is in most cases not possible due to the unknown errors in the input parameters. By studying the sensitivity of the commonly and widely used conventional BHT correction methods of Horner plot, Lachenbruch and Brewer method, Brennand’s method, forward modeling, linearization method, and 1BHT method, we developed a workflow to adjust BHT to SFT considering parameter availability and uncertainty. Usually, the BHT-corrected values are given as defined temperatures, ignoring the error that lies within the input parameters. Instead, our method aims to provide probability scenarios (e.g., p10, p50, and p90 limits) that can be used as expected value, worst-case, and best-case scenario and that can be used as a business case for the successful implementation of geothermal projects. Since the thermal output of a hydrothermal well depends on production temperatures, the large uncertainty ranges that can occur are a serious concern when estimating the expectable risk to economic efficiency. In addition, the probability scenarios can be used in the evaluation of borehole data, e.g., in the correction of hydraulic parameters or calculation of hydro-geothermal parameters for which the borehole temperature must be known.

In the future, a regional representation of predicted formation temperatures can build on this approach as it allows for a valid propagation of the likely errors to be considered. Furthermore, the error-based correction workflow can be transferred to other geothermal settings when tested on known static formation temperatures that can be estimated from drill stem tests, fiber-optic temperature sensing, or wireline temperature logs.