- Email: [email protected]
Effects of topography and climatic changes on the temperature in borehole GFU-1, Prague
TECTONOPHYSICS ELSEVIER
Tectonophysics 239 (1994) 187-197
Effects of topography and climatic changes on the temperature in borehole GFU-1, Prague Jan Safanda Geophysical Institute of the Academy of Sciences of the Czech Republic, Bohi II, 14131 Prague, Czech Republic
Received 11 April 1994; revised version accepted
18 July 1994
Abstract In October, 1992, a 150 m borehole was drilled on the ground of the Geophysical Institute, Prague, for purposes of long-term temperature observations. Since the borehole is intended for a detailed study of the change during downward propagation of the ground surface temperature (GST), an evaluation of the topographic effect and the unsteady-state component is necessary. As a first step, the effect of topography was studied by solving the heat conduction equation in 2-D models. The GST model presumed a lapse rate of 5 mK/m, superposed on the GST dependence on the slope orientation, which raises the GST on the southern slope by 0.5 K and lowers it by the same amount on the northern one. As a second step, the model was subjected to 1 K warming of the surface 70 yr ago (which is in the range of expected changes) and the computed profile was inverted by a 1-D algorithm to reconstruct
the GST. Neglecting the influence of topography, the reconstructed warming gives a value of 1.05 K 74 yr ago. In view of these results, the inversion of the measured profile was made without applying the topographic correction. If a step change in the GST is presumed, a warming of 1.2 K 70 yr ago is most probable; an assumption of a linear change gives a warming of 1.1 K 120 yr ago. The meteorological data recorded at Praha-Klementinum support the latter case.
1. Introduction Since the very beginning of heat flow measurements it has been known that subsurface temperature gradients, determined primarily by the amount of heat coming to the surface from the deeper parts of the Earth, are also influenced by effects of shallow origin. From the point of view of the thermal state of the Earth’s interior, the influence of such factors as topography, groundwater convection, sedimentation, erosion, longterm ground surface temperature (GST) changes, etc., produce disturbances which must be removed. The last factor, the GST changes, has 0040-195 l/94/$07.00 0 1994 Elsevier Science SSDI 0040.1951(94)00146-4
been known for many years. One of the first to recognize the significant effect of the past glaciation history on the subsurface temperature distribution was Lane (1923), followed by Benfield (19391, Bullard (19391, Birch (1948) and many others. Gradually, it was realized that the corrections should be applied not only to the areas of past glaciation, but worldwide and that younger climatic events should also be taken into consideration (Beck, 1977; Zoth and Haenel, 1988). The solution of the opposite problem of inferring past climatic variations from subsurface temperature data was started, with the exception of an isolated approach of Hotchkiss and Ingersoll
B.V. All rights reserved
J. &fanda / Tectonophysics 239 (1994) 187-197
188
Kubik, 1992). One of the crucial questions is the coupling between climatic and GST changes. In order to understand this process more clearly, a 150 m deep borehole, referred to as GFU-1, was drilled in October, 1992, in the grounds of the Geophysical Institute in Prague. Another bore-
(19341, by Beck and Judge (1969), followed by Cermak (1971), Anderssen and Saul1 (1973) and Beck (19821, among others. At present, the problem is being studied from many aspects (Beck, 1992; Beck et al., 1992; Glow, 1992j Lewis and Wang, 1992; Shen and Beck, 1992; Safanda and 0
PARALLEL
NORMAL
++
a)
++
0
a
+
1
l
+
j++ I +
+++
+
+ ++
t
+ +
.‘* .’
.
l
. +
’ .
+
l
*
.
. I
l
.
l
.*
~
‘*’
++ ++ + ++ ,/I .+* : ++ !-_-AL +
.
+
+/‘++
++
1 -1
150
+
+
++++ +++ +
+
.
(+ \ ++
+
5. l. .
l
l.
*
+++
1
100
+
f ++ + ++*
1
96
+* f *.****. + +++ +
+
+
*.
*
+ # ++* +
50-
a
++
+’
. + + + + ++
‘.
+
l
1.0
1.5
2.b
2.b
3.b
conductivity,
3.b
(4
W/mK
0 NORMAL
VERTICAL
PARALLEIL
b\
___
150
1.0
1.5
2.0
2.5 conductivity,
Fig. 1. Thermal conductivity in borehole GFU-1. (a) Components parallel 3rd-order polynomials. (b) Vertical component and the best fit approximation
3.0
3.5
4
W/mK and normal to the bedding. The best fit lines are of all three components by 7th-order polynomials.
J. .&f&da / Ttctonophysics 239 (1994) 187-197
hole was drilled a few metres apart to the depth of 40 m. The shallower hole was equipped with a set of thermistors for permanent monitoring of the temperature. The data will be combined with repeated logging of the GFU-1 hole and with ground and air temperature measurements carried out in the same regime as at a meteorological station. The observation programme should provide a sufficient data set for analyzing the problem in its full complexity. The paper deals with the first stage of the project, which consisted in interpreting the equilibrium profile of GFU-1.
2. Data set Borehole GFU-1 is situated in the oldest and most stable part of the Bohemian Massif, which belongs to the Variscan branch of the Hercynian system. It encountered sediments of Ordovician age, belonging to the Zahoiany (O-126 m) and Vinice (126-150 m) Members, mostly represented by micaceous, silty to clayey shales. The drill core recovery was almost lOO%, with the exception of the depth interval 38-48 m, where some parts of the core were broken to small chips. The diameter of the core was 108 mm. The maximum inclination of the borehole, 2.1”, was measured at a depth of 100 m. The inclination of the bedding varied between 0” and 90” as a result of folding. The mean value was about 50” to the north. The thermal conductivity, both along the bedding and perpendicular (normal) to it, was measured on core samples from 124 depths by the commercial Quick Thermal Meter (QTM-2) (Sass et al., 1983; CermLk et al., 1984; Galson et al., 1987) based on a transient line heat source method. The samples were cut perpendicular to the bedding. In such a case, both the above-mentioned principal components of the conductivity tensor can be determined from two measurements, during which the QTM instrument lies on the cut sample, first parallel to the bedding and then perpendicular to it (Grubbe et al., 1983). All samples were processed and their conductivity measured l-3 days after the core recovery. During the measurement, the samples were placed in
189
a vessel filled with water-saturated sand, so that only 1-2 cm of them were protruding. Their surface was covered by a 0.01 mm waterproof foil in order to prevent evaporation. The conductivity values obtained in this way are believed to be very close to the in situ ones (Fig. la). The mean value from all 124 depths is 3.2 k 0.2 W/(mK) along the bedding and 2.2 * 0.4 W/(mK) normal to it. This means that the anisotropy factor is 1.45. Besides the two principal conductivities measured, the vertical conductivity was computed at each of 124 depths using the formula for the transformation of conductivity tensor depending on the inclination of the bedding (Carslaw and Jaeger, 1959) (Fig. lb). The vertical conductivity, which determines in a decisive manner the vertical temperature gradient, varies by 3 W/(mK) in the upper 100 m and displays a pronounced minimum of about 2.2 W/(mK) at depths of 120-130 m. The minimum is the consequence of the decrease in the inclination of the bedding, combined with the decrease in both principal conductivities. The radioactive heat production was measured on 10 core samples, more or less regularly spaced throughout the whole borehole. Samples of about 2 kg in weight were crushed ( < 0.5 mm) and used for the gamma-spectrometric measurement of Th, U and K contents. The analysis was carried out with the 8192-channel gamma-spectrometer Cicero (Silena, Milano) (Vaiikova et al., 1993). Values of the Th, U and K concentrations obtained were used for calculating the heat production by means of the equation given in KreSl and Vaiikova (1978). The production ranges between 1.70 X 10e6 and 2.16 X lo-” W/m” and the mean value is (2.0 + 0.1) X 10Ph W/m-‘. The temperature in the borehole was measured by a portable thermometer (KreSl, 1981). The readings were taken point by point at 5 m intervals. The accuracy of the probe can be estimated to be better than f0.05 K and the relative precision of the temperature readings is better than 0.01 K. An estimate of the equilibrium temperature was determined by analysis of 11 loggings made 2-78 days after drilling (Stulc, 1994). The mean annual surface temperature extrapolated from the depth interval of 20-30 m is 9.9”C. This is 0.6 K higher than the value resulting from
J. hfanda
190
/ Tectonophysics 239 (1994) 187-197
encountered, 2~10~~ W/m3, accounts for the difference of 0.3 mW/m2 between the surface and the bottom of the hole only.
the regression formula (Kubik, 1990), based on data from 48 hydrometeorological stations in the Czech Republic. The vertical temperature gradient computed by means of the central difference operator for intervals of 10 m is shown in Fig. 2. The values oscillate appreciably but the general trend is an increase from about lo-15 mK/m in the first 50 m to 20-25 mK/m in the lowermost 50 m. In explaining the observed increase in the temperature gradient with depth, we can rule out two main causes; namely the inclination of the borehole and groundwater movement. The inclination is negligible and there are no indications of groundwater activity. The picture does not change after computing the heat flow in the 10 m intervals (Fig. 2). The mean vertical conductivity values seem to be representative enough since they are based on 8-11 measurements in each interval. Despite this, the heat flow oscillates between 31 and 43 mW/m2 in the uppermost 50 m and between 50 and 63 mW/m* below 100 m. It should be noticed that the heat production conductivity,
1
GRADIENT
3. Effect of topography Another possible reason for heat flow increasing with depth is the effect of topography. The borehole is located on a flat elevation extending for a few hundred metres in an E-W direction. The flank toward the south is steep, about ll”, and only 20 m high, whereas the flank to the north drops 60 m at an angle of < 5”. In a first approximation the structure can be treated as a 2-D one. Thus, the effect of topography was studied by solving numerically the heat conduction equation in 2-D geothermal models of the borehole surroundings (Safanda, 1987). In order to take into account the strong anisotropy of the thermal conductivity, the computer program described in Safanda (1994) was used. The GST W/mK
CONDUCTIVITY HEAT
150 ~““,““,‘“‘,‘“‘,“,,,,,“,,~,,,,,,,,,,,,,”’,,,,,,,,,,i 10 15 20 25 30 35 gradient,
K/km;
FLOW
_ 40
45
50
55
60
_ 65
i
heat flow, mW/m2
Fig. 2. Temperature gradient, vertical conductivity and heat flow in GFU-1 for the 10 m intervals The heavy lines of the gradient and heat flow represent values corrected for the topography.
in the lo-150
m depth
section.
191
J. .fafanda / Tectonophysics 239 (1994) 187-197 5120
4S20
4720
4520
4320
4120
~9~,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,~~,,,,,,,,,,,,,,,,,,,~~~~~~~
0
100
150
300 3920
4320
4120
4720
4520
--DISTANCE,
5120
4920
300
M
Fig. 3. Effect of topography at the site of GFU-1 shown in the form of isolines of the vertical heat flow normalized to the basal heat flow at a depth of 5.5 km. The southern flank is 0.5 K warmer and the northern flank is 0.5 K cooler than the flat terrain. An anisotropic conductivity structure inclines 50” to the north.
normalized 0.
0
E
1
0.40
0.50
0.60
vertical 0.70
heat 0.80
flow
-8 0.90
1.00
50
1.1
- 50
4 5
a
100
- 100
150
0.30
0.40
0.50
0.60
0.70
0.80
0.90
1.00
1.10
150
Fig. 4. Distortion of the heat flow in borehole GFU-1 due to the topography for different models of the GST and thermal conductivity. la = the results of the preferred model, also shown in Fig. 3, the GST depends on the slope, anisotropic layers inclined 50” to the north; lb = as la but layers are inclined 50” to the south; 2a = the GST does not depend on the slope, an anisotropic conductivity model, layers inclined to the north; 2b = as 2a but layers inclined to the south; 3 = the GST does not depend on the slope, an isotropic conductivity model; 4 = the GST depends on the slope, an isotropic conductivity model; 5 = the GST depends on the slope, an isotropic conductivity model, the GST increased by 1 K 70 yr ago.
192
J. hfanda
/ Tectonophysics
model presumed an atmospheric lapse rate of 5 mK/m (Kubik, 1990). According to Blackwell et al. (19801, the GST is dependent, at least in the middle latitudes, on the angle and orientation of the slope. The results published applied to our topography yield a difference in GST of 1.2-1.5 K between the southern and northern slopes at the same altitude. That is why the northern flank was kept 0.5 K cooler and the southern 0.5 K warmer than would be expected from the lapse rate considered. As it is evident that the pattern of folding is local, due to its small scale, the conductivity model was based on a mean bedding inclination of 50 to the north and the mean conductivity 3.2 W/(mK) along the bedding and 2.2 W/(mK) perpendicular to it. The boundary condition at the bottom, at a depth of 5.5 km, was a uniform heat flow of 60 mW/m2, which is believed to be close to the undisturbed value. At the sides of the model, 4.5 km away from the borehole, conditions of horizontal symmetry of the temperature field were prescribed. The heat sources within the model were put equal to zero in order not to obscure the pure effect of the topography. Results of the modelling are shown in Fig. 3 in a form of isolines of the vertical heat flow normalized to the undisturbed basal value. Since the depth is reckoned from the highest point of the model surface (outside the figure), the head of the well has a coordinate of 30 m. The heat flow pattern near the surface is rather complicated, which is given mainly by the schematic approximation of the surface. Nevertheless, values of the normalized heat flow at the site of the borehole increase steadily from 0.86 at 25 m, through 0.90 at 100 m, to 0.93 at 150 m. A region of the disturbance less than 0.5% starts at the depth of about 500 m. The results say that the measured heat flow at 25 m is 14% lower than the undisturbed one, at 100 m 10% lower and at 150 m 7%. Some other models were considered, too. The results are summarized in Fig. 4. Apart from the model described above (curve la), a model with the same anisotropic rock but with the opposite inclination of the strata, 50” to the south (curve lb), was studied. It shows a distortion up to 40%
239 (1994) 187-197
lower than curve la. The difference can be ascribed to the asymmetric surface temperature on the northern and southern flanks adjacent to the plateau carrying the borehole. This conclusion is based on the results displayed by curves 2a and 2b. They represent the same conductivity models as la and lb, respectively, but the surface temperature depends on the altitude only. For the isotropic conductivity, the dependence of the heat flow distortion on the surface temperature model is shown by curves 3 (the GST is a function of the altitude only) and 4 (the GST depends on the orientation of the slope). It can be concluded that the most realistic model, la, yields the greatest heat flow distortion in the depth interval used for the terrestrial heat flow determination (below 20 m). The model was used to correct the observed values of the vertical gradient and heat flow. As shown in Fig. 2, the corrected curves are approximately parallel with the uncorrected ones; that is the correction does not remove the trend of their increase with depth. The magnitude of the shift is between +5 and +8 mW/m2. The mean value of the corrected heat flow in the loo-150 m section is 62 L- 6 mW/m2.
4. Effect of the change in ground surface temperature
Since it was found that the topography cannot account for the heat flow increase with depth, we investigated the influence of a possible long-term change in the GST. As mentioned in the introduction, this topic is being studied intensively and borehole GFU-1 is intended for such research. The algorithm used for the inversion of the observed temperature profile in terms of the GST change was described in Safanda and Kubik (1992). To evaluate the character of the change, heat conduction in a homogeneous halfspace was presumed. Only step and linear changes in the GST were considered. The parameters estimated by the inversion are: the magnitude of the GST change and the time of its inception, together with the mean annual GST and the subsurface gradient prior to the change. The advantage of the method is the possibility of quantifying our
J. iafanda
/ Tectonophysics
in the a priori estimates of the parameters sought in a form of their a priori standard deviations (SDS) and the possibility of obtaining estimates of the a posteriori SDS. The a posteriori to a priori SDS ratio (SDR) is used to characterize the reduction in the uncertainty of the parameters estimated. As a first step, the method was applied to the measured profile. Because the temperature, and not the temperature gradient, is inverted, the knowledge of the gradient topographic correction cannot be used to separate the effect of the topography and the change in GST. In processing the profile, the O-25 m section, influenced by annual variation, was omitted. The values resulting from the linear regression applied to the lowermost section below 125 m were taken as an a priori estimate of the undisturbed temperature gradient and surface temperature. A priori SDS of the magnitude and time of the GST change were chosen broad enough to leave the parameters effectively unconstrained and easily adjusted by the inversion to satisfy the measured temperature, which is much more strongly constrained. The magnitude of the a priori SD of one temper-
Fig. 5. The temperature 1.2 K 70 yr ago.
by
profile
from borehole
GFU-1
193
239 (1994) 187-197
ature measurement, 0.02 K, was twice as much as the estimated experimental error. The a posteriori parameters of a step change are a warming by 1.2 K (SDR = 0.3) 70 yr ago (SDR = 0.41, with an undisturbed GST of 8.8”C (SDR = 0.3) and a gradient of 24.5 mK/m (SDR = 0.3). If a linear change in GST is assumed, the inversion indicated a warming by 1.1 K (SDR = 0.2) 120 yr ago (SDR = 0.8), with an undisturbed GST of 9.o”C (SDR = 0.2) and a gradient of 23.3 mK/m (SDR = 0.3). The results show that the most poorly determined quantity is the time of onset of the change, which displays a SDR between 0.4 and 0.8. While the difference in the magnitude of the change in the GST between the step and linear models is lo%, the times differ by more than 40%. The mean-square difference between the measured and computed temperature curves is 0.017 K for the step change and 0.021 K for the linear one. The degree of coincidence is shown in graphic form in Fig. 5. The theoretical curves for the two types of GST change considered merge on the scale of the figure. The question now is, how are the results bi-
and the theoretical
curve computed
as a response
to the increase
in the GST
194
J. &fmda
/ Tectonophysics 239 (1994) 187-197
ased by interpreting the measured temperature curve entirely as a response to the GST change without a correction for the effect of topography. To assess the share of the two effects, the warming was simulated mathematically by solving the unsteady-state heat conduction equation in the 2-D geothermal model of the borehole surroundings, which was identical with that used in calculating the topographic corrections. The steadystate temperature distribution corresponding to the topographic correction model was used as an initial condition. The temperature-depth profiles calculated are shown in Fig. 6. They represent the response of the initial steady-state profile to 1 K warming at various times after the onset of the change. In contrast to the real data, the simulation allows the unsteady-state component of the temperature to be extracted from the ‘measured’ curve. For instance, the difference between curve 70 and curve ‘steady-state’ in Fig. 6 represents the unsteady-state component induced by 1 K step warming 70 yr ago. Its inversion yields a warming of 1.05 K 74 yr ago. This small deviation
from the true values arises from differences in the change propagation in the 2-D geothermal model used in the temperature calculation and the 1-D theory used in the inversion algorithm. To simulate the real situation, when the measured profile cannot be deprived of the topographic effect, synthetic curve 70 from Fig. 6 was inverted without subtracting the steady-state component. Results for the step change indicate warming by 1.13 K (SDR = 0.3) 66 yr ago (SDR = 0.3). Thus, the estimated GST increase is 13% higher and the time of the change 6% lower than the true values. In this case the discrepancy issues not only from the differences between 2-D and 1-D models used for the direct and inverse problems, respectively, but also from an inclusion of the steady-state topographic effect into the inversion. If the steady-state curve from Fig. 6, alone, was inverted, the estimated parameters of the step change suggested a warming by 0.13 K (SDR = 0.1) 44 yr ago (SDR = 0.8). These results confirm the impression obtained from visual inspection of Fig. 6, where the curvature of the steady-
Fig. 6. Temperature-depth profiles in GFU-1 obtained by solving the heat conduction equation in the 2-D model of the borehole surroundings. The steady-state curve reflects the pure effect of the topography. The other curves represent the response of the steady-state curve to a 1 K step increase in the GST 10, 30, 50 and 70 yr after the change and for infinite time.
195
J. .?afanda / Tectonophysics 239 (1994) 187-197
state curve, appreciably smaller than that of the unsteady-state curves, suggests a secondary role for the surrounding topography in the evaluation of the change in the GST. The magnitude of the warming, 0.13 K, exactly accounts for the surplus of the warming observed in the previous inversion, where a ‘sum’ of the unsteady- and steadystate components was interpreted. Nevertheless, the exact coincidence of the two values is accidental because the inversion of the pure unsteady-state component does not yield 1 K warming but 1.05 K, as presented above. It can be concluded that the analysis of the measured, uncorrected GFU-1 temperature profile in terms of the change in GST slightly overestimates (by O.l0.2 K) the magnitude of the warming and negligibly underestimates (by a few years) the time at which it begins.
5. Conclusions An analysis of the equilibrium temperature profile in borehole GFU-1, together with detailed thermal conductivity measurements, has shown that the increase in the temperature gradient with depth is not compensated for by decreasing conductivity. In addition, the effect of the topography, evaluated by solving the heat conduction equation in a 2-D geothermal model of the borehole surroundings, can only explain a smaller part of the increase. It was found that its decisive part can be ascribed to the GST warming by about 1 K. The time of the change depends on the shape considered and varies between 70 yr ago for a step change and 120 yr ago for a linear increase. The unsteady-state solution of the heat conduction equation, in the 2-D model simulating the simultaneous effect of the topography and the GST warming, can be used as a correction of the heat flow observed. The heat flow in the borehole, normalised to the value undisturbed by the topography and the GST change, attains 0.364 for a 1 K warming 70 yr ago at a depth of 25 m, 0.735 at 100 m and 0.890 at 150 m. The mean value, calculated as an average from corrected values at the 10 m intervals, is 82 + 17 mW/m2. The mean heat flow, corrected for the topogra-
6.5
”
10.0
-
PRAHA-KLEMENTINUM
6.0
/
E
firl
E : 9.5 k
5.5
3
5.0
9.0 -
6.5”“““““““““““’ 1750 1600
1650 years,
1900 A.D.
1950
Fig. 7. The 10 yr averages of the air temperature Prague-Klementinum (the left-hand temperature at Milesovka (the right-hand scale).
4.5 2000
measured at scale) and
phy only, amounts to 62 + 6 mW/m’ in the 70150 m section. The increase due to the climatic correction represents more than 30%. Since borehole GFU-1 is expected to be open for temperature logging for many years, it will be possible to monitor the time changes of the subsurface temperature field caused by the GST increase propagating downward. The theoretical maximum annual change induced by the warming by 1 K 70 yr ago, 0.0035 K, should appear at a depth of 60-70 m. Such a temperature difference cannot be detected reliably by the present equipment. However, a cumulative effect could enlarge the difference to a detectable value within a few years. In addition, a new temperature probe with greater accuracy is being constructed, which promises to record even the predicted annual changes. Information about the warming can be compared with direct instrumental measurements of the air temperature which have been recorded at Prague-Klementinum for more than 220 yr (Fig. 7). The records (PetroviE, 1969; Jilek, 1990) reveal a permanent warming trend since the middle of the 19th century. The amplitude of the change, calculated as the difference between the mean temperature of the warmest decade (1980-1990), +10.0.X and the coldest decade (1850-18601, +8.7”C, is + 1.35 K over a timespan of 130 yr.
196
J. Shfanda / Tectonophysics 239 (1994) 187-197
These characteristics are very close to those obtained by the GFU-1 profile inversion, 1.1 K 120 yr ago, if a linear shape to the change in the GST is considered. It should be noted, however, that the borehole is located on the rim of Prague, about 7 km from the centre, where the Klementinum station is situated. Part of the warming observed at Klementinum could be caused by the growth of the city, which accelerated in the second half of last century. The 10 yr averages of the air temperature measured at Milesovka station since 1905 are shown in Fig. 7 for comparison. The Milesovka station is located 60 km northwest of Prague at an altitude of 800 m. Unlike Prague, it is out of the reach of any direct anthropogenic influence and its records may indicate a regional warming trend. The warming does not seem to proceed so rapidly at this station. The increase during the period 1905-1990 here is 0.65 K instead of the 0.95 K recorded at Klementinum. It is to be hoped that at least some of the questions discussed will be answered by detailed monitoring of the borehole.
Acknowledgements This study was performed within projects 205/93/0412 and 205/93/0413 funded by the Grant Agency of the Czech Republic. The author would like to express his thanks to Dr. Vera VankovB, who provided data on the radioactive heat production and to Mgr. Petr Stulc, who made most of the temperature loggings.
References Anderssen, R.S. and Saull, V.A., 1973. Surface temperature history determinations from borehole measurements. Math. Geol., 5: 269-283. Beck, A.E., 1977. Climatically perturbed temperature gradients and their effect on regional and continental heat flow means. Tectonophysics, 41: 17-39. Beck, A.E., 1982. Precision logging of temperature gradients and the extraction of past climate. Tectonophysics, 83: l-11. Beck, A.E., 1992. Inferring past climate change from subsur-
face profiles: some problems and methods. Global Planet. Change, 98: 73-80. Beck, A.E. and Judge, A.S., 1969. Analyses of heat flow data: detailed observations in a single borehole. Geophys. J. R. Astron. Sot., 18: 145-158. Beck, A.E., Shen, P.Y., Beltrami, H., Mareshal, J.-C., Safanda, J., Sebagenzi, M.N., Vasseur, G. and Wang, K., 1992. A comparison of five different analyses in the interpretation of five data sets. Global Planet. Change, 98: 101-112. Benfield, A.E., 1939. Terrestrial heat flow in Great Britain. Proc. R. Sot. London Ser. A, 173: 428-450. Birch, F., 1948. The effects of Pleistocene climatic variations upon geothermal gradients. Am. J. Sci., 246: 729-760. Blackwell, D.D., Steele, J.L. and Brott, C.A., 1980. The terrain effect on terrestrial heat flow. J. Geophys. Res., 85: 4757-4772. Bullard, E.C., 1939. Heat flow in South Africa. Proc. R. Sot. London Ser. A, 173: 474-502. Carslaw, H.S. and Jaeger, J.C., 1959. Conduction of Heat in Solids. Clarendon Press, Oxford, 510 pp. Cermik, V., 1971. Underground temperatures and inferred climatic temperature of the past millenium. Palaeogeogr. Palaeoclimatol. Palaeoecol., 10: 1-19. Cermik, V., KreSl, M., Safanda, J., Napoles-Pruna, M., Tenreyro-Perez, R., Torres, L.M. and Valdes, J.J., 1984. First heat flow density assessment in Cuba. Tectonophysics, 103: 283-296. Glow, G.D., 1992. The extent of temporal smearing in surface-temperature histories derived from borehole temperature measurements. Global Planet. Change, 98: 81-86. Galson, D.A., Wilson, N.P., Scharli, U. and Rybach, L., 1987. A comparison of the divided-bar and QTM methods of measuring thermal conductivity. Geothermics, 16(3): 215226. Grubbe, K., Haenel, R. and Zoth, G., 1983. Determination of the vertical component of thermal conductivity by line source methods. Zentralbl. Geol. Paliontol. Teil 1, l/2: 49-56. Hotchkiss, W.O. and Ingersoll, L.R., 1934. Post-glacial time calculations from recent geothermal measurements in the Calumet Copper Mines. J. Gee\, 42: 113-142. Jilek, J., 1990. Statisticka roEenka Ceski a Slovenskt Federativni Republiky 1990. SNTL, Prague, 740 pp. Krell, M., 1981. Technique of heat flow measurements. In: A. Zatopek (Editor), Geophysical Synthesis in Czechoslovakia. Veda, Bratislava, pp. 449-454. KreSI, M. and Vafikova, V., 1978. Radioactivity and heat production data from several boreholes in the Bohemian Massif. Stud. Geophys. Geod., 22: 165-176. Kubik, J., 1990. Subsurface temperature field of the Bohemian Massif. Stud. Geophys. Geod., 34: 110-128. Lane, AC., 1923. Geotherms of Lake Superior Copper Country. Bull. Geol. Sot. Am., 34: 703-720. Lewis, T.J. and Wang, K., 1992. Influence of terrain on bedrock temperatures. Global Planet. Change, 48: 87-100. PetroviE, S., 1969. Podnebi Ceskoslovenske socialisticke republiky. Hydrometeorol. Inst., Prague, 356 pp.
J. Shfanda / Tectonophysics 239 (I 994) 187-l 97 Safanda, J., 1987. Some remarks on the estimation of geothermal topocorrections. Stud. Geophys. Geod., 31: 284-300. Safanda, J., 1994. Effect of thermal conductivity anisotropy of rocks on the subsurface temperature field. Geophys. J. Int. (submitted). Safanda, J. and Kubik, .I., 1992. Evidence of ground surface temperature changes from two boreholes in the Bohemian Massif. Global Planet. Change, 98: 199-208. Sass, J.H., Stone, C. and Munroe, R.J., 1983. Thermal conductivity determination on solid rock-a comparison between the steady-state divided bar apparatus and a commercial transient line-source device. J. Volcanol. Geotherm. Res., 20: 145-153. Shen, P.Y. and Beck, A.E., 1992. Paleoclimate change and
197
heat flow density inferred from temperature data in the Superior Province of the Canadian Shield. Global Planet. Change, 98: 143-165. Stulc, P., 1994. Temperature recovery after drilling in the borehole GFU-1 in Prague. Tectonophysics (submitted) Vaiikova, V., Vanek, J., HanuS, V., Fiala, J. and Chalupslj, V., 1993. Geochemistry of talc-alkaline volcanic rocks in relation to deep structure of Ecuador and Southern Colombia. Global Tectonics Metallogeny, 4(3): 113-155. Zoth, G. and Haenel, R., 1988. Paleotemperatures. In: R. Haenel, L. Rybach and L. Stegena (Editors), Handbook of Terrestrial Heat-Flow Density Determinations. Kluwer, Dordrecht, Appendix 10.2.
Tectonophysics 239 (1994) 187-197
Effects of topography and climatic changes on the temperature in borehole GFU-1, Prague Jan Safanda Geophysical Institute of the Academy of Sciences of the Czech Republic, Bohi II, 14131 Prague, Czech Republic
Received 11 April 1994; revised version accepted
18 July 1994
Abstract In October, 1992, a 150 m borehole was drilled on the ground of the Geophysical Institute, Prague, for purposes of long-term temperature observations. Since the borehole is intended for a detailed study of the change during downward propagation of the ground surface temperature (GST), an evaluation of the topographic effect and the unsteady-state component is necessary. As a first step, the effect of topography was studied by solving the heat conduction equation in 2-D models. The GST model presumed a lapse rate of 5 mK/m, superposed on the GST dependence on the slope orientation, which raises the GST on the southern slope by 0.5 K and lowers it by the same amount on the northern one. As a second step, the model was subjected to 1 K warming of the surface 70 yr ago (which is in the range of expected changes) and the computed profile was inverted by a 1-D algorithm to reconstruct
the GST. Neglecting the influence of topography, the reconstructed warming gives a value of 1.05 K 74 yr ago. In view of these results, the inversion of the measured profile was made without applying the topographic correction. If a step change in the GST is presumed, a warming of 1.2 K 70 yr ago is most probable; an assumption of a linear change gives a warming of 1.1 K 120 yr ago. The meteorological data recorded at Praha-Klementinum support the latter case.
1. Introduction Since the very beginning of heat flow measurements it has been known that subsurface temperature gradients, determined primarily by the amount of heat coming to the surface from the deeper parts of the Earth, are also influenced by effects of shallow origin. From the point of view of the thermal state of the Earth’s interior, the influence of such factors as topography, groundwater convection, sedimentation, erosion, longterm ground surface temperature (GST) changes, etc., produce disturbances which must be removed. The last factor, the GST changes, has 0040-195 l/94/$07.00 0 1994 Elsevier Science SSDI 0040.1951(94)00146-4
been known for many years. One of the first to recognize the significant effect of the past glaciation history on the subsurface temperature distribution was Lane (1923), followed by Benfield (19391, Bullard (19391, Birch (1948) and many others. Gradually, it was realized that the corrections should be applied not only to the areas of past glaciation, but worldwide and that younger climatic events should also be taken into consideration (Beck, 1977; Zoth and Haenel, 1988). The solution of the opposite problem of inferring past climatic variations from subsurface temperature data was started, with the exception of an isolated approach of Hotchkiss and Ingersoll
B.V. All rights reserved
J. &fanda / Tectonophysics 239 (1994) 187-197
188
Kubik, 1992). One of the crucial questions is the coupling between climatic and GST changes. In order to understand this process more clearly, a 150 m deep borehole, referred to as GFU-1, was drilled in October, 1992, in the grounds of the Geophysical Institute in Prague. Another bore-
(19341, by Beck and Judge (1969), followed by Cermak (1971), Anderssen and Saul1 (1973) and Beck (19821, among others. At present, the problem is being studied from many aspects (Beck, 1992; Beck et al., 1992; Glow, 1992j Lewis and Wang, 1992; Shen and Beck, 1992; Safanda and 0
PARALLEL
NORMAL
++
a)
++
0
a
+
1
l
+
j++ I +
+++
+
+ ++
t
+ +
.‘* .’
.
l
. +
’ .
+
l
*
.
. I
l
.
l
.*
~
‘*’
++ ++ + ++ ,/I .+* : ++ !-_-AL +
.
+
+/‘++
++
1 -1
150
+
+
++++ +++ +
+
.
(+ \ ++
+
5. l. .
l
l.
*
+++
1
100
+
f ++ + ++*
1
96
+* f *.****. + +++ +
+
+
*.
*
+ # ++* +
50-
a
++
+’
. + + + + ++
‘.
+
l
1.0
1.5
2.b
2.b
3.b
conductivity,
3.b
(4
W/mK
0 NORMAL
VERTICAL
PARALLEIL
b\
___
150
1.0
1.5
2.0
2.5 conductivity,
Fig. 1. Thermal conductivity in borehole GFU-1. (a) Components parallel 3rd-order polynomials. (b) Vertical component and the best fit approximation
3.0
3.5
4
W/mK and normal to the bedding. The best fit lines are of all three components by 7th-order polynomials.
J. .&f&da / Ttctonophysics 239 (1994) 187-197
hole was drilled a few metres apart to the depth of 40 m. The shallower hole was equipped with a set of thermistors for permanent monitoring of the temperature. The data will be combined with repeated logging of the GFU-1 hole and with ground and air temperature measurements carried out in the same regime as at a meteorological station. The observation programme should provide a sufficient data set for analyzing the problem in its full complexity. The paper deals with the first stage of the project, which consisted in interpreting the equilibrium profile of GFU-1.
2. Data set Borehole GFU-1 is situated in the oldest and most stable part of the Bohemian Massif, which belongs to the Variscan branch of the Hercynian system. It encountered sediments of Ordovician age, belonging to the Zahoiany (O-126 m) and Vinice (126-150 m) Members, mostly represented by micaceous, silty to clayey shales. The drill core recovery was almost lOO%, with the exception of the depth interval 38-48 m, where some parts of the core were broken to small chips. The diameter of the core was 108 mm. The maximum inclination of the borehole, 2.1”, was measured at a depth of 100 m. The inclination of the bedding varied between 0” and 90” as a result of folding. The mean value was about 50” to the north. The thermal conductivity, both along the bedding and perpendicular (normal) to it, was measured on core samples from 124 depths by the commercial Quick Thermal Meter (QTM-2) (Sass et al., 1983; CermLk et al., 1984; Galson et al., 1987) based on a transient line heat source method. The samples were cut perpendicular to the bedding. In such a case, both the above-mentioned principal components of the conductivity tensor can be determined from two measurements, during which the QTM instrument lies on the cut sample, first parallel to the bedding and then perpendicular to it (Grubbe et al., 1983). All samples were processed and their conductivity measured l-3 days after the core recovery. During the measurement, the samples were placed in
189
a vessel filled with water-saturated sand, so that only 1-2 cm of them were protruding. Their surface was covered by a 0.01 mm waterproof foil in order to prevent evaporation. The conductivity values obtained in this way are believed to be very close to the in situ ones (Fig. la). The mean value from all 124 depths is 3.2 k 0.2 W/(mK) along the bedding and 2.2 * 0.4 W/(mK) normal to it. This means that the anisotropy factor is 1.45. Besides the two principal conductivities measured, the vertical conductivity was computed at each of 124 depths using the formula for the transformation of conductivity tensor depending on the inclination of the bedding (Carslaw and Jaeger, 1959) (Fig. lb). The vertical conductivity, which determines in a decisive manner the vertical temperature gradient, varies by 3 W/(mK) in the upper 100 m and displays a pronounced minimum of about 2.2 W/(mK) at depths of 120-130 m. The minimum is the consequence of the decrease in the inclination of the bedding, combined with the decrease in both principal conductivities. The radioactive heat production was measured on 10 core samples, more or less regularly spaced throughout the whole borehole. Samples of about 2 kg in weight were crushed ( < 0.5 mm) and used for the gamma-spectrometric measurement of Th, U and K contents. The analysis was carried out with the 8192-channel gamma-spectrometer Cicero (Silena, Milano) (Vaiikova et al., 1993). Values of the Th, U and K concentrations obtained were used for calculating the heat production by means of the equation given in KreSl and Vaiikova (1978). The production ranges between 1.70 X 10e6 and 2.16 X lo-” W/m” and the mean value is (2.0 + 0.1) X 10Ph W/m-‘. The temperature in the borehole was measured by a portable thermometer (KreSl, 1981). The readings were taken point by point at 5 m intervals. The accuracy of the probe can be estimated to be better than f0.05 K and the relative precision of the temperature readings is better than 0.01 K. An estimate of the equilibrium temperature was determined by analysis of 11 loggings made 2-78 days after drilling (Stulc, 1994). The mean annual surface temperature extrapolated from the depth interval of 20-30 m is 9.9”C. This is 0.6 K higher than the value resulting from
J. hfanda
190
/ Tectonophysics 239 (1994) 187-197
encountered, 2~10~~ W/m3, accounts for the difference of 0.3 mW/m2 between the surface and the bottom of the hole only.
the regression formula (Kubik, 1990), based on data from 48 hydrometeorological stations in the Czech Republic. The vertical temperature gradient computed by means of the central difference operator for intervals of 10 m is shown in Fig. 2. The values oscillate appreciably but the general trend is an increase from about lo-15 mK/m in the first 50 m to 20-25 mK/m in the lowermost 50 m. In explaining the observed increase in the temperature gradient with depth, we can rule out two main causes; namely the inclination of the borehole and groundwater movement. The inclination is negligible and there are no indications of groundwater activity. The picture does not change after computing the heat flow in the 10 m intervals (Fig. 2). The mean vertical conductivity values seem to be representative enough since they are based on 8-11 measurements in each interval. Despite this, the heat flow oscillates between 31 and 43 mW/m2 in the uppermost 50 m and between 50 and 63 mW/m* below 100 m. It should be noticed that the heat production conductivity,
1
GRADIENT
3. Effect of topography Another possible reason for heat flow increasing with depth is the effect of topography. The borehole is located on a flat elevation extending for a few hundred metres in an E-W direction. The flank toward the south is steep, about ll”, and only 20 m high, whereas the flank to the north drops 60 m at an angle of < 5”. In a first approximation the structure can be treated as a 2-D one. Thus, the effect of topography was studied by solving numerically the heat conduction equation in 2-D geothermal models of the borehole surroundings (Safanda, 1987). In order to take into account the strong anisotropy of the thermal conductivity, the computer program described in Safanda (1994) was used. The GST W/mK
CONDUCTIVITY HEAT
150 ~““,““,‘“‘,‘“‘,“,,,,,“,,~,,,,,,,,,,,,,”’,,,,,,,,,,i 10 15 20 25 30 35 gradient,
K/km;
FLOW
_ 40
45
50
55
60
_ 65
i
heat flow, mW/m2
Fig. 2. Temperature gradient, vertical conductivity and heat flow in GFU-1 for the 10 m intervals The heavy lines of the gradient and heat flow represent values corrected for the topography.
in the lo-150
m depth
section.
191
J. .fafanda / Tectonophysics 239 (1994) 187-197 5120
4S20
4720
4520
4320
4120
~9~,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,~~,,,,,,,,,,,,,,,,,,,~~~~~~~
0
100
150
300 3920
4320
4120
4720
4520
--DISTANCE,
5120
4920
300
M
Fig. 3. Effect of topography at the site of GFU-1 shown in the form of isolines of the vertical heat flow normalized to the basal heat flow at a depth of 5.5 km. The southern flank is 0.5 K warmer and the northern flank is 0.5 K cooler than the flat terrain. An anisotropic conductivity structure inclines 50” to the north.
normalized 0.
0
E
1
0.40
0.50
0.60
vertical 0.70
heat 0.80
flow
-8 0.90
1.00
50
1.1
- 50
4 5
a
100
- 100
150
0.30
0.40
0.50
0.60
0.70
0.80
0.90
1.00
1.10
150
Fig. 4. Distortion of the heat flow in borehole GFU-1 due to the topography for different models of the GST and thermal conductivity. la = the results of the preferred model, also shown in Fig. 3, the GST depends on the slope, anisotropic layers inclined 50” to the north; lb = as la but layers are inclined 50” to the south; 2a = the GST does not depend on the slope, an anisotropic conductivity model, layers inclined to the north; 2b = as 2a but layers inclined to the south; 3 = the GST does not depend on the slope, an isotropic conductivity model; 4 = the GST depends on the slope, an isotropic conductivity model; 5 = the GST depends on the slope, an isotropic conductivity model, the GST increased by 1 K 70 yr ago.
192
J. hfanda
/ Tectonophysics
model presumed an atmospheric lapse rate of 5 mK/m (Kubik, 1990). According to Blackwell et al. (19801, the GST is dependent, at least in the middle latitudes, on the angle and orientation of the slope. The results published applied to our topography yield a difference in GST of 1.2-1.5 K between the southern and northern slopes at the same altitude. That is why the northern flank was kept 0.5 K cooler and the southern 0.5 K warmer than would be expected from the lapse rate considered. As it is evident that the pattern of folding is local, due to its small scale, the conductivity model was based on a mean bedding inclination of 50 to the north and the mean conductivity 3.2 W/(mK) along the bedding and 2.2 W/(mK) perpendicular to it. The boundary condition at the bottom, at a depth of 5.5 km, was a uniform heat flow of 60 mW/m2, which is believed to be close to the undisturbed value. At the sides of the model, 4.5 km away from the borehole, conditions of horizontal symmetry of the temperature field were prescribed. The heat sources within the model were put equal to zero in order not to obscure the pure effect of the topography. Results of the modelling are shown in Fig. 3 in a form of isolines of the vertical heat flow normalized to the undisturbed basal value. Since the depth is reckoned from the highest point of the model surface (outside the figure), the head of the well has a coordinate of 30 m. The heat flow pattern near the surface is rather complicated, which is given mainly by the schematic approximation of the surface. Nevertheless, values of the normalized heat flow at the site of the borehole increase steadily from 0.86 at 25 m, through 0.90 at 100 m, to 0.93 at 150 m. A region of the disturbance less than 0.5% starts at the depth of about 500 m. The results say that the measured heat flow at 25 m is 14% lower than the undisturbed one, at 100 m 10% lower and at 150 m 7%. Some other models were considered, too. The results are summarized in Fig. 4. Apart from the model described above (curve la), a model with the same anisotropic rock but with the opposite inclination of the strata, 50” to the south (curve lb), was studied. It shows a distortion up to 40%
239 (1994) 187-197
lower than curve la. The difference can be ascribed to the asymmetric surface temperature on the northern and southern flanks adjacent to the plateau carrying the borehole. This conclusion is based on the results displayed by curves 2a and 2b. They represent the same conductivity models as la and lb, respectively, but the surface temperature depends on the altitude only. For the isotropic conductivity, the dependence of the heat flow distortion on the surface temperature model is shown by curves 3 (the GST is a function of the altitude only) and 4 (the GST depends on the orientation of the slope). It can be concluded that the most realistic model, la, yields the greatest heat flow distortion in the depth interval used for the terrestrial heat flow determination (below 20 m). The model was used to correct the observed values of the vertical gradient and heat flow. As shown in Fig. 2, the corrected curves are approximately parallel with the uncorrected ones; that is the correction does not remove the trend of their increase with depth. The magnitude of the shift is between +5 and +8 mW/m2. The mean value of the corrected heat flow in the loo-150 m section is 62 L- 6 mW/m2.
4. Effect of the change in ground surface temperature
Since it was found that the topography cannot account for the heat flow increase with depth, we investigated the influence of a possible long-term change in the GST. As mentioned in the introduction, this topic is being studied intensively and borehole GFU-1 is intended for such research. The algorithm used for the inversion of the observed temperature profile in terms of the GST change was described in Safanda and Kubik (1992). To evaluate the character of the change, heat conduction in a homogeneous halfspace was presumed. Only step and linear changes in the GST were considered. The parameters estimated by the inversion are: the magnitude of the GST change and the time of its inception, together with the mean annual GST and the subsurface gradient prior to the change. The advantage of the method is the possibility of quantifying our
J. iafanda
/ Tectonophysics
in the a priori estimates of the parameters sought in a form of their a priori standard deviations (SDS) and the possibility of obtaining estimates of the a posteriori SDS. The a posteriori to a priori SDS ratio (SDR) is used to characterize the reduction in the uncertainty of the parameters estimated. As a first step, the method was applied to the measured profile. Because the temperature, and not the temperature gradient, is inverted, the knowledge of the gradient topographic correction cannot be used to separate the effect of the topography and the change in GST. In processing the profile, the O-25 m section, influenced by annual variation, was omitted. The values resulting from the linear regression applied to the lowermost section below 125 m were taken as an a priori estimate of the undisturbed temperature gradient and surface temperature. A priori SDS of the magnitude and time of the GST change were chosen broad enough to leave the parameters effectively unconstrained and easily adjusted by the inversion to satisfy the measured temperature, which is much more strongly constrained. The magnitude of the a priori SD of one temper-
Fig. 5. The temperature 1.2 K 70 yr ago.
by
profile
from borehole
GFU-1
193
239 (1994) 187-197
ature measurement, 0.02 K, was twice as much as the estimated experimental error. The a posteriori parameters of a step change are a warming by 1.2 K (SDR = 0.3) 70 yr ago (SDR = 0.41, with an undisturbed GST of 8.8”C (SDR = 0.3) and a gradient of 24.5 mK/m (SDR = 0.3). If a linear change in GST is assumed, the inversion indicated a warming by 1.1 K (SDR = 0.2) 120 yr ago (SDR = 0.8), with an undisturbed GST of 9.o”C (SDR = 0.2) and a gradient of 23.3 mK/m (SDR = 0.3). The results show that the most poorly determined quantity is the time of onset of the change, which displays a SDR between 0.4 and 0.8. While the difference in the magnitude of the change in the GST between the step and linear models is lo%, the times differ by more than 40%. The mean-square difference between the measured and computed temperature curves is 0.017 K for the step change and 0.021 K for the linear one. The degree of coincidence is shown in graphic form in Fig. 5. The theoretical curves for the two types of GST change considered merge on the scale of the figure. The question now is, how are the results bi-
and the theoretical
curve computed
as a response
to the increase
in the GST
194
J. &fmda
/ Tectonophysics 239 (1994) 187-197
ased by interpreting the measured temperature curve entirely as a response to the GST change without a correction for the effect of topography. To assess the share of the two effects, the warming was simulated mathematically by solving the unsteady-state heat conduction equation in the 2-D geothermal model of the borehole surroundings, which was identical with that used in calculating the topographic corrections. The steadystate temperature distribution corresponding to the topographic correction model was used as an initial condition. The temperature-depth profiles calculated are shown in Fig. 6. They represent the response of the initial steady-state profile to 1 K warming at various times after the onset of the change. In contrast to the real data, the simulation allows the unsteady-state component of the temperature to be extracted from the ‘measured’ curve. For instance, the difference between curve 70 and curve ‘steady-state’ in Fig. 6 represents the unsteady-state component induced by 1 K step warming 70 yr ago. Its inversion yields a warming of 1.05 K 74 yr ago. This small deviation
from the true values arises from differences in the change propagation in the 2-D geothermal model used in the temperature calculation and the 1-D theory used in the inversion algorithm. To simulate the real situation, when the measured profile cannot be deprived of the topographic effect, synthetic curve 70 from Fig. 6 was inverted without subtracting the steady-state component. Results for the step change indicate warming by 1.13 K (SDR = 0.3) 66 yr ago (SDR = 0.3). Thus, the estimated GST increase is 13% higher and the time of the change 6% lower than the true values. In this case the discrepancy issues not only from the differences between 2-D and 1-D models used for the direct and inverse problems, respectively, but also from an inclusion of the steady-state topographic effect into the inversion. If the steady-state curve from Fig. 6, alone, was inverted, the estimated parameters of the step change suggested a warming by 0.13 K (SDR = 0.1) 44 yr ago (SDR = 0.8). These results confirm the impression obtained from visual inspection of Fig. 6, where the curvature of the steady-
Fig. 6. Temperature-depth profiles in GFU-1 obtained by solving the heat conduction equation in the 2-D model of the borehole surroundings. The steady-state curve reflects the pure effect of the topography. The other curves represent the response of the steady-state curve to a 1 K step increase in the GST 10, 30, 50 and 70 yr after the change and for infinite time.
195
J. .?afanda / Tectonophysics 239 (1994) 187-197
state curve, appreciably smaller than that of the unsteady-state curves, suggests a secondary role for the surrounding topography in the evaluation of the change in the GST. The magnitude of the warming, 0.13 K, exactly accounts for the surplus of the warming observed in the previous inversion, where a ‘sum’ of the unsteady- and steadystate components was interpreted. Nevertheless, the exact coincidence of the two values is accidental because the inversion of the pure unsteady-state component does not yield 1 K warming but 1.05 K, as presented above. It can be concluded that the analysis of the measured, uncorrected GFU-1 temperature profile in terms of the change in GST slightly overestimates (by O.l0.2 K) the magnitude of the warming and negligibly underestimates (by a few years) the time at which it begins.
5. Conclusions An analysis of the equilibrium temperature profile in borehole GFU-1, together with detailed thermal conductivity measurements, has shown that the increase in the temperature gradient with depth is not compensated for by decreasing conductivity. In addition, the effect of the topography, evaluated by solving the heat conduction equation in a 2-D geothermal model of the borehole surroundings, can only explain a smaller part of the increase. It was found that its decisive part can be ascribed to the GST warming by about 1 K. The time of the change depends on the shape considered and varies between 70 yr ago for a step change and 120 yr ago for a linear increase. The unsteady-state solution of the heat conduction equation, in the 2-D model simulating the simultaneous effect of the topography and the GST warming, can be used as a correction of the heat flow observed. The heat flow in the borehole, normalised to the value undisturbed by the topography and the GST change, attains 0.364 for a 1 K warming 70 yr ago at a depth of 25 m, 0.735 at 100 m and 0.890 at 150 m. The mean value, calculated as an average from corrected values at the 10 m intervals, is 82 + 17 mW/m2. The mean heat flow, corrected for the topogra-
6.5
”
10.0
-
PRAHA-KLEMENTINUM
6.0
/
E
firl
E : 9.5 k
5.5
3
5.0
9.0 -
6.5”“““““““““““’ 1750 1600
1650 years,
1900 A.D.
1950
Fig. 7. The 10 yr averages of the air temperature Prague-Klementinum (the left-hand temperature at Milesovka (the right-hand scale).
4.5 2000
measured at scale) and
phy only, amounts to 62 + 6 mW/m’ in the 70150 m section. The increase due to the climatic correction represents more than 30%. Since borehole GFU-1 is expected to be open for temperature logging for many years, it will be possible to monitor the time changes of the subsurface temperature field caused by the GST increase propagating downward. The theoretical maximum annual change induced by the warming by 1 K 70 yr ago, 0.0035 K, should appear at a depth of 60-70 m. Such a temperature difference cannot be detected reliably by the present equipment. However, a cumulative effect could enlarge the difference to a detectable value within a few years. In addition, a new temperature probe with greater accuracy is being constructed, which promises to record even the predicted annual changes. Information about the warming can be compared with direct instrumental measurements of the air temperature which have been recorded at Prague-Klementinum for more than 220 yr (Fig. 7). The records (PetroviE, 1969; Jilek, 1990) reveal a permanent warming trend since the middle of the 19th century. The amplitude of the change, calculated as the difference between the mean temperature of the warmest decade (1980-1990), +10.0.X and the coldest decade (1850-18601, +8.7”C, is + 1.35 K over a timespan of 130 yr.
196
J. Shfanda / Tectonophysics 239 (1994) 187-197
These characteristics are very close to those obtained by the GFU-1 profile inversion, 1.1 K 120 yr ago, if a linear shape to the change in the GST is considered. It should be noted, however, that the borehole is located on the rim of Prague, about 7 km from the centre, where the Klementinum station is situated. Part of the warming observed at Klementinum could be caused by the growth of the city, which accelerated in the second half of last century. The 10 yr averages of the air temperature measured at Milesovka station since 1905 are shown in Fig. 7 for comparison. The Milesovka station is located 60 km northwest of Prague at an altitude of 800 m. Unlike Prague, it is out of the reach of any direct anthropogenic influence and its records may indicate a regional warming trend. The warming does not seem to proceed so rapidly at this station. The increase during the period 1905-1990 here is 0.65 K instead of the 0.95 K recorded at Klementinum. It is to be hoped that at least some of the questions discussed will be answered by detailed monitoring of the borehole.
Acknowledgements This study was performed within projects 205/93/0412 and 205/93/0413 funded by the Grant Agency of the Czech Republic. The author would like to express his thanks to Dr. Vera VankovB, who provided data on the radioactive heat production and to Mgr. Petr Stulc, who made most of the temperature loggings.
References Anderssen, R.S. and Saull, V.A., 1973. Surface temperature history determinations from borehole measurements. Math. Geol., 5: 269-283. Beck, A.E., 1977. Climatically perturbed temperature gradients and their effect on regional and continental heat flow means. Tectonophysics, 41: 17-39. Beck, A.E., 1982. Precision logging of temperature gradients and the extraction of past climate. Tectonophysics, 83: l-11. Beck, A.E., 1992. Inferring past climate change from subsur-
face profiles: some problems and methods. Global Planet. Change, 98: 73-80. Beck, A.E. and Judge, A.S., 1969. Analyses of heat flow data: detailed observations in a single borehole. Geophys. J. R. Astron. Sot., 18: 145-158. Beck, A.E., Shen, P.Y., Beltrami, H., Mareshal, J.-C., Safanda, J., Sebagenzi, M.N., Vasseur, G. and Wang, K., 1992. A comparison of five different analyses in the interpretation of five data sets. Global Planet. Change, 98: 101-112. Benfield, A.E., 1939. Terrestrial heat flow in Great Britain. Proc. R. Sot. London Ser. A, 173: 428-450. Birch, F., 1948. The effects of Pleistocene climatic variations upon geothermal gradients. Am. J. Sci., 246: 729-760. Blackwell, D.D., Steele, J.L. and Brott, C.A., 1980. The terrain effect on terrestrial heat flow. J. Geophys. Res., 85: 4757-4772. Bullard, E.C., 1939. Heat flow in South Africa. Proc. R. Sot. London Ser. A, 173: 474-502. Carslaw, H.S. and Jaeger, J.C., 1959. Conduction of Heat in Solids. Clarendon Press, Oxford, 510 pp. Cermik, V., 1971. Underground temperatures and inferred climatic temperature of the past millenium. Palaeogeogr. Palaeoclimatol. Palaeoecol., 10: 1-19. Cermik, V., KreSl, M., Safanda, J., Napoles-Pruna, M., Tenreyro-Perez, R., Torres, L.M. and Valdes, J.J., 1984. First heat flow density assessment in Cuba. Tectonophysics, 103: 283-296. Glow, G.D., 1992. The extent of temporal smearing in surface-temperature histories derived from borehole temperature measurements. Global Planet. Change, 98: 81-86. Galson, D.A., Wilson, N.P., Scharli, U. and Rybach, L., 1987. A comparison of the divided-bar and QTM methods of measuring thermal conductivity. Geothermics, 16(3): 215226. Grubbe, K., Haenel, R. and Zoth, G., 1983. Determination of the vertical component of thermal conductivity by line source methods. Zentralbl. Geol. Paliontol. Teil 1, l/2: 49-56. Hotchkiss, W.O. and Ingersoll, L.R., 1934. Post-glacial time calculations from recent geothermal measurements in the Calumet Copper Mines. J. Gee\, 42: 113-142. Jilek, J., 1990. Statisticka roEenka Ceski a Slovenskt Federativni Republiky 1990. SNTL, Prague, 740 pp. Krell, M., 1981. Technique of heat flow measurements. In: A. Zatopek (Editor), Geophysical Synthesis in Czechoslovakia. Veda, Bratislava, pp. 449-454. KreSI, M. and Vafikova, V., 1978. Radioactivity and heat production data from several boreholes in the Bohemian Massif. Stud. Geophys. Geod., 22: 165-176. Kubik, J., 1990. Subsurface temperature field of the Bohemian Massif. Stud. Geophys. Geod., 34: 110-128. Lane, AC., 1923. Geotherms of Lake Superior Copper Country. Bull. Geol. Sot. Am., 34: 703-720. Lewis, T.J. and Wang, K., 1992. Influence of terrain on bedrock temperatures. Global Planet. Change, 48: 87-100. PetroviE, S., 1969. Podnebi Ceskoslovenske socialisticke republiky. Hydrometeorol. Inst., Prague, 356 pp.
J. Shfanda / Tectonophysics 239 (I 994) 187-l 97 Safanda, J., 1987. Some remarks on the estimation of geothermal topocorrections. Stud. Geophys. Geod., 31: 284-300. Safanda, J., 1994. Effect of thermal conductivity anisotropy of rocks on the subsurface temperature field. Geophys. J. Int. (submitted). Safanda, J. and Kubik, .I., 1992. Evidence of ground surface temperature changes from two boreholes in the Bohemian Massif. Global Planet. Change, 98: 199-208. Sass, J.H., Stone, C. and Munroe, R.J., 1983. Thermal conductivity determination on solid rock-a comparison between the steady-state divided bar apparatus and a commercial transient line-source device. J. Volcanol. Geotherm. Res., 20: 145-153. Shen, P.Y. and Beck, A.E., 1992. Paleoclimate change and
197
heat flow density inferred from temperature data in the Superior Province of the Canadian Shield. Global Planet. Change, 98: 143-165. Stulc, P., 1994. Temperature recovery after drilling in the borehole GFU-1 in Prague. Tectonophysics (submitted) Vaiikova, V., Vanek, J., HanuS, V., Fiala, J. and Chalupslj, V., 1993. Geochemistry of talc-alkaline volcanic rocks in relation to deep structure of Ecuador and Southern Colombia. Global Tectonics Metallogeny, 4(3): 113-155. Zoth, G. and Haenel, R., 1988. Paleotemperatures. In: R. Haenel, L. Rybach and L. Stegena (Editors), Handbook of Terrestrial Heat-Flow Density Determinations. Kluwer, Dordrecht, Appendix 10.2.