Determination of Seafloor Temperatures Using Data from High-Resolution Marine Heat Probes

1 Institute of Oceanography, National Taiwan University, Taipei, Taiwan, ROC * Corresponding author address: Prof. Chuen-Tien Shyu, Institute of Oceanography, National Taiwan University, Taipei, Taiwan, ROC; E-mail: ctshyu@oc.ntu.edu.tw Marine heat flow surveys typically utilize Lister-type heat probes, which allow multi-penetrations to measure both temperature and thermal conductivity in situ. The temperature rise resulting from penetration friction must first be removed in order to obtain the equilibrium temperature and the thermal conductivity of the sediments. However, the time needed to return to the original sediment temperature may be very long. Using highresolution heat probes, temperatures with resolution of 0.1 mK were obtained within the first few meters below the seafloor off southwestern Taiwan. This research introduces an improved data processing technique that can extrapolate friction-raised temperature to infinite time and obtain equilibrium ambient temperature accurately from a limited amount of highresolution temperature recordings. The computing technique is to regress the temperature decay function numerically on parameters of ambient temperature, initial temperature rise, thermal conductivity, and heat capacity of the sediment. In order to obtain a better fit between the data and the model as well as a better prediction of the friction raised temperature, we allowed erroneous solutions to parameters for regression. This technique was tested on both synthetic and high-resolution field data to demonstrate accuracy and efficiency. From the field experiment, we estimated that the error between the prediction and the true equilibrium ambient temperatures to be within a few 0.1 mK which is much smaller than a typical resolution, of a few mK in marine heat flow measurement. (


INTRODUCTION
Heat flow on the seafloor is determined by independent measurements of the vertical temperature gradient and thermal conductivity using a heat probe in near surface sediments.In order to make these measurements, high-resolution sensors must be utilized to resolve small differences between any two sensors spaced vertically about 50 to 100 cm.The penetration of the probe introduces transient disturbances to the temperature of the surrounding sediments through friction.In principle, it takes a significant amount of time for the temperature to return to the original undisturbed sediment temperature.For example, we need to keep the heat probe undisturbed in the bottom sediment for about four hours if we expect to obtain the equilibrium temperature of the sediment accurate to 1 mK (10 3 − K) and for forty hours to 0.1 mK (10 4 − K) according to Eq.1 (assuming friction-induced temperature rise T 0 = 0.3°C, heat capacity S = 2.5 10 Jm K 6 × − − 3 1 , radius a = 4 mm, and thermal conductivity K = 1 Wm K − − 1 1 ).Due to disturbances caused either by coring or ambient environmental conditions (e.g., pressure, temperature, water content), conductivity measurements in situ appear to compare favorably with those made on core material either aboard ship or in the laboratory (Davis et al. 1984).Von Herzen et al.(1982) and Jemsek et al. (1985) developed the in situ method with outrigged probes, using a line heater located within the same probe as the thermistors used for gradient determinations.When continuous heating was applied for specified time periods (e.g., 10 minutes), a significant fraction of the applied heat had dissipated to the sediment after the recording period.The disadvantage of using too much power becomes particularly evident when more measurements are desired at a site for each heat probe that is lowered from the vessel.Lister (1979) proposed a pulse heating method as an alternative to the continuous heating process.In his method, a short calibrated heat pulse is applied to a cylindrical probe, and analytic approximations are used to fit the decay of the heat pulse in order to determine conductivity (Lister 1979;Hyndman et al. 1979).A temperature-time record of a typical penetration is depicted in Fig. 1.This record indicates three operational stages: (1) the time period prior to the final descent of the probe and penetration, during which constant depth (i.e., constant temperature) is maintained (isothermal stage), (2) the thermal decay which resulted from the impact of the probe on the sediments (impact decay) and (3) the thermal decay following the imposition of an artificial heat pulse for conductivity measurement (heat-pulse decay).The relative offset of each temperature sensor was determined during the isothermal stage.The temperature and temperature gradient were determined by extrapolating the temperature of the impact decay (about 12 minutes) to infinite time.However, extrapolation may be unreliable beyond the recorded time range.Therefore, the extrapolation error must be accounted for properly.
The temperature rising from penetration friction affects the original ambient temperature of the sediments as well as the successive heat pulse experiment, thereby affecting the measurement of temperature gradients and thermal conductivities.Although the notion of removing the effect of frictional heating from the field data through infinite extrapolation has been mentioned in the literature, the procedure has seldom been explicitly presented.Based on the decay function for a cylindrical probe, Davis et al. (1984) described an iterative data reduction scheme to calculate heat flow at four sites on the Blake-Bahama outer ridge.They concluded that thermal conductivities measured in situ provide the most reliable heat flow determinations.Villinger and Davis (1987) discussed a similar data reduction algorithm in detail.It was concluded that, according to the numerical experiment, the algorithm they developed may be able to determine the undisturbed sediment temperature to within 1 mK and the in situ thermal conductivity may be accurate to within 1%.Shyu (2000) presented an algorithm for extrapolating the penetration decay temperature to infinite time to obtain the undisturbed sediment temperature.The decay curve was fitted with the cylindrical decay function through regression with respect to ambient temperature, friction-induced temperature rise, and the diffusivity of the sediment.Based on a similar approach, Hartmann and Villinger (2002) have proposed a different calculation step for data processing; they have achieved a standard deviation of 1.4 mK (0.28 per cent) for ambient temperature and 0.018 Wm K − − 1 1 (1.8 per cent) for thermal conductivity.
Due to instrumental and ambient noises, it is commonly believed that the temperature uncertainties are slightly greater than 1 mK when using a heat probe to measure heat flows on the ocean floor (Ge'li et al. 2001;Fisher et al. 2002;Hartmann and Villinger 2002).However, Fig. 1.A typical example of temperature vs. time during a heat flow measurement using a pulse heating probe with 7 thermistors.
a more accurate or better resolution is needed in some specific cases such as detecting the base of gas hydrate stability zones.If high resolution can be achieved successfully, it is advantageous to use the pulse-heating type probe, which requires substantially less energy, takes less heating duration and allows a smaller sampling interval.The effect of shortening the heating duration is equivalent to compressing the heat pulse into a spike-like wavelet, which is very useful for the post-data processing.For example, it makes the determination of the temperature rise time easier.On the other hand, in order to make the instrument construction easier and the heat flow measurements on the ocean floor more efficient, the idea of considering the frictional heat to be a heat pulse for determining thermal properties has been proposed for years by Lee and Von Herzen (1994).Recently, Lee et al. (2003) further presented an inverse modeling method, which formulated together with the finite-element simulation of the friction-heat cooling history, and concluded that the equilibrium temperature can be determined within 0.022 K and the determination of conductivity and diffusivity typically around 5 per cent.However, one of the inherent drawbacks of implementing this idea is that the amplitude of frictional temperature-rise is very small (in general, less than 0.5°C mainly depending on the length of the probe and sediment properties); this may be partly overcome using high-resolution data.As a whole, whatever the improvement of the probe construction, the data reduction technique described above would make the heat flow measurements on the ocean floor more efficient, save ship-operation costs, and make data processing easier.
Based on a heat-flow probe capable of resolving 0.1 mK in the range of -1°C to 25°C, we present numerical experiments and field examples to demonstrate the use of such a highresolution marine heat probe.We propose a computational technique to obtain a better fit between the high-resolution data and the model, as well as a better prediction of equilibrium ambient temperature and the extrapolation temperature from data of short frictional temperature recordings.

REDUCTION ALGORITHM FOR SUCCESSIVE TEMPERATURE EXTRAPOLATION
Assuming that a perfect conducting cylindrical probe with an initial temperature T 0 is inserted into the sediment with an undisturbed temperature T B lower than T 0 and that there is no contact resistance, the temperature of the probe T at a time t is (Bullard 1954;Jaeger 1956;Carslaw and Jaeger 1959) where J u n ( ) and Y u n ( ) are respectively Bessel functions of order n of the first and second kinds.τ = kt a / 2 is the dimensionless time constant with the probe of radius a in a material of diffusivity k, k K c = / ρ is the ratio of thermal conductivity to heat capacity of the sediment, α π ρ = 2 2 a c/S is twice the ratio of the heat capacity ( π ρ c a 2 ) per unit length of the sediment to that of the probe, ρ and c are the density and specific heat of the sediment, respectively, and S is the heat capacity of the probe.
Equation 1 is nonlinear with respect to K and ρc (linear with respect to T B and T 0 ).Our algorithm is based on the Gauss-Newton method which minimizes the sum of squares of the misfits between data and expectations.The key concept underlying the technique is a Taylor series expansion, which is used to express the original nonlinear equation in an approximate, linear form.Then, the least-squares method is employed to obtain new estimates of the parameters that move toward minimizing the misfits.Since S and a of the probe are known, after expanding the nonlinear Eq. 1 in a Taylor series around parameters, T B , T 0 , K and ρc the first derivates can be expressed in the matrix form: where vector D { } contains the differences between the measurements and the Eq. 1 values at

M , and Z j
[ ] is the matrix of partial derivatives of the Eq. 1 at the jth iteration where n is the number of data points.The vector ∆V { } contains the changes in the param- eter values The column vector E { } contains the residuals e e e e n 1 2 3 L { } T and: The partial derivatives of the Eq. 1 with respect to K and ρc are difficult to evaluate.We use difference equations to approximate the partial derivatives as: where h K and h ρc are specified small increments of K and ρc, respectively.According to the linear least-squares theory, the sum of the squares of the residuals can be minimized by taking its partial derivative with respect to each of the coefficients and setting the equation equal to zero.The outcome of this process is the normal equations which can be expressed concisely in matrix form as: Thus, the approach consists of solving Eq. 3 for ∆V { }, which is employed iteratively to improve the parameter as in V V V i j i j i , , This procedure is repeated until the solution converges-that is, until However, if we solve the four parameters simultaneously from Eq. 3, three possible shortcomings may emerge.( 1 Jm K (Hyndman et al. 1979), and then to calculate T 0 , K, and ρc with the previous solved T B at the second stage.The equilibrium ambient tempera- ture T B is an independent variable in Eq. 1, which can be accurately solved during the first computing stage.

PENETRATION
The sensor tube used in the study was constructed mainly of steel, mineral oil, heater wires, Teflon insulated cables, and thermistors.Since these materials are not "perfect conductors", there is a time delay in sensor response to frictional heating.The initial temperature rise for each thermistor is different (approximately from 0.1 to 0.5°C), depending on the sediment properties and the penetration depth.The variation in effective origin times for each thermistor is substantial due to changing properties and varied sediment disturbances along the length of the probe and the sampling intervals.Actually, frictional heating is time-varied, but the initial temperature T 0 in Eq. 1 is assumed to be a spike-like impulse.For the application, we may contract the initially raised friction temperature T 0 into a spike with a delay time.Moreover, the variation in effective origin times for the frictional dissipation is generally small (within 10 seconds).During the computation, the iteration loops are run with different delays (e.g., multiple of sampling interval).Because the misfits are very sensitive to the delays, the effective origin time can be determined with an accuracy of < 10 seconds (i.e., two sampling intervals) readily if the data sampling interval is 5 seconds.Through the time delay correction, the remaining errors resulting from time-shift may be absorbed mainly by T 0 , K, and ρc during regression and will not affect the temperature extrapolation, which is discussed in the following examples.

SYNTHETIC DATA
The first example is to examine the deviation of the temperature gradient caused by frictional heating.Assuming the surface temperature on the seafloor is T surface = 4 °C, the ambient temperature gradient of the sediment is G = 0.0350°C m −1 , the thermal conductivity of the sediment is K = 0.9 Wm K − − 1 1 , and the heat capacity of the sediment is 3 31 10 6 3 1 . × − − Jm K , after the impact of the probe, the temperature rises for sensor 1 and sensor 2 are T 01 = 0.3°C and T 02 = 0.15°C, respectively.The interval between the two sensors is 0.7 m.The simulated temperature gradient and the given value are shown in Fig. 2. It is noted that their discrepancy is larger than 10% at 12 minutes.We observe that this discrepancy curve asymptotically approaches zero at infinite time.
The second illustration which, tests the accuracy of extrapolated temperature from synthetic data of different resolutions, examines a heat probe impacted on the sediment with an ambient temperature T B = 4.0°C and thermal conductivity K = 1.0 Wm K − − 1 1 .The friction raises temperature to T 0 = 0.3°C after the impact, and the temperatures after penetration were recorded for 12 minutes.The temporal sampling interval is 5 seconds.We consider synthetic temperatures to have different resolutions of 0.00001, 0.1, and 1 mK, respectively.That is, digits after 0.00001, 0.1, and 1 mK of the synthetic data are omitted.In addition, the effective origin time for the recording may be uncertain by as much as 10 seconds.The synthetic data from 7 to 12 minutes with different temperature resolutions and time shifts are fitted with Eq. 1 by regressing against variables T B , T 0 and K first and then T 0 , K, and ρc as the proposed method; the calculated results are given in Table 1.
After performing the regression, we can compute Sr, the sum of the squares of the residuals around the regression line and St, the total sum of the squares around the mean for the dependent variable (in our case, the temperatures).The r (St -Sr)/ St 2 = is called the coefficient of determination and r is the correlation coefficient.In all test cases, the high values correlation coefficients (≅ 1) indicate that the fits are very good.However, the resulting T 0 Table 1.Comparison the calculated parameters with those of the given model.a The resolution of the temperature data.b Correlation coefficient.and K were deviated from the input values used to derive the synthetic data.This implies that it is not feasible to resolve the temperature decay function into variables T 0 , K, and ρc despite the high resolving power in the estimate of ambient equilibrium temperatures, T B .Though, it is unrealistic to expect a measured temperature on the ocean floor to achieve a resolution level of 10 8 − K, the solutions for these variables obtained from five-digit resolution (0.1 mK, highresolution) data does exhibit a great increase in accuracy over that from four-digit resolution (1 mK) data.The extrapolation from a short duration (7 to 12 minutes) of high-resolution data is consistent with the synthetic data (0.00001 mK resolution) that lie beyond 12 minutes (Fig. 3).The ambient equilibrium temperature and the extrapolated temperature obtained by the method highly depend upon data resolution.The results for delayed origin time (+10 seconds) are relatively better than that of the advanced origin time (-10 seconds) (Fig. 4).Nevertheless, the maximum error for extrapolations, up to 2 hours for both cases, is less than 0.1 m°C.This is better than the level of instrumental resolution.The results strongly indicate that the extrapolating temperature and the equilibrium ambient temperature T B of the sediment can be solved simultaneously with high accuracy from the transient decay data, even if the effective origin time and the magnitude of friction raised temperature are uncertain.The algorithm described above is particularly useful for the determination of the ambient equilibrium temperature and the extrapolated temperature of the sediment from a short duration of high-resolution frictional temperature decay data.It differs from the iterative scheme proposed by Hartmann and Villinger (2002) in that it employs T B , T 0 , K, and ρc to regress the temperature decay function with different time delays, while Hartmann and Villinger (2002) employs frictional heat, equilibrium temperature (T B ), conductivity (K), and time delay to regress the temperature decay function.Their examples on synthetic data show that equilibrium temperature can be determined with an error of about 1~2 mK, or with a relative error of 0.28% which is more than one order of magnitude greater than that of the high-resolution data and, therefore, is not suitable for processing high-resolution data.

EXAMPLES OF FIELD DATA
In situ heat flow measurements were made at a station southwest off Taiwan.During the experiment, a Lister-type probe developed at the Institute of Oceanography, National Taiwan University (Fig. 5) was used.It allows temperature resolution to be better than 0.1 mK in the range of -1 to 25°C.Concerning the effectiveness of the high-resolution probe, we performed two penetrations with different resolution probes at the same location .Figure 6 shows one example demonstrating the temperature decays measured at about 3 meters below the seafloor (mbsf) (water depth 1460 m) resulting from frictional heating after the entry of probes into sediment.A distinct difference is observed between the temperatures obtained with a typical heat probe with temperature resolution of 1 mK (upper panel in Fig. 6) and those obtained with 0.1 mK (lower panel in Fig. 6).Clearly, perturbations arising from oceanographic disturbances (e.g., bottom water temperature fluctuations, geothermal ambient noise) at the test location are not significant (< 0.1 mK) at depth of 3~4 mbsf (lower panel in Fig. 6).We note that the large deviations from the regression curve (upper panel in Fig. 6) and call into question the accuracy realized with typical probes with resolution of 1 mK commonly used in heat-flow studies.
Testing the goodness of the extrapolation from a short duration of high-resolution data were conducted; Fig. 7 shows a typical data set.When the sensor string penetrates the sedi- Inside the tube, thermistor spacing is about 60~100 cm (depending on the length of the tube).Total weight is about 350 kg.ment the friction heat raising the temperature decay was recorded at a sample rate of 5 sec for 2.4 hours (2840 to 11480 sec) until the probe was pulled out of the sediment.Seven minutes of data duration from 5 to 12 minutes after penetration (at 2840 sec) were used for the extrapolation (up to 2.88 hours after penetration).The standard error of the estimate was 6 4 10 5 .× − °C after 12 iterations for thermistor sensor 1 and 6 3 10 5 .× − °C after 11 iterations for thermistor sensor 2, respectively.These errors are acceptable since the resolution of the high-resolution marine heat probe used for the measurement is about 0.1 mK.The corresponding equilibrium ambient temperatures 2.7770°C and 2.7345°C for each sensor were solved simultaneously.The extrapolated temperature beyond 12 minutes after penetration can be used to remove the residual frictional heat effect during the following onset of heat pulse decay used for thermal conductivity reduction.The maximum deviation of about 0.1 m°C between the extrapolated temperature and that of the measured data occurred at 11400 seconds or about 2.38 hours after the penetration which conformed to the requirement.

CONCLUTIONS
A pulse-heating probe was widely used to measure heat flow on the ocean floor in order to save energy.There were two periods of temperature response (impact decay and heat-pulse decay) after penetration by the thermistors constructed inside the tube of the probe.In order to minimize the effect of penetration heat disturbance when obtaining temperature, temperature gradient, and thermal conductivity of the sediment, it was necessary to extend the measuring time significantly until the penetration heat dissipated.This is not always feasible or possible, especially when considering constraints such as cruise time and total operation cost.One practical aspect is to predict accurately the friction-raised temperatures to infinite time by using an efficient extrapolation computing technique performed on a short period of high-resolution data.
Contrary to the common belief that the ambient temperature noise on the ocean floor is about 1 mK, we found that the ambient temperature in the first few mbsf at a site off southwestern Taiwan with a water depth of 1460 m to be relatively stable (temperature noise < 0.1 m°C) and was suitable for high-resolution temperature (e.g., 0.1 m°C) measurements.Furthermore, in order to meet the requirements for the proposed data reduction technique, the heat probe should provide reading resolution up to four decimal places in °C.
Since high resolution measurements are allowed to sample sufficiently finely, so that the uncertainty of the effective origin time can be continually adjusted during the regression process until it has converged within 10 seconds.The errors caused by time-shift may be absorbed by T 0 , K, and ρc (especially by K and ρc) during regression and will not affect the tempera- ture extrapolation (accurate to few 0.1 m°C).
The temperature decay that resulted from the penetration frictional heat may be adequately described with an instantaneously heated cylindrical probe.By fitting the temperature decay curve (about 12 minute duration) with Eq. 1, and by regressing with respect to parameters T B , T 0 and K at first calculation stage, one can obtain accurate equilibrium ambient temperature T B as well as accurate temperature gradients.Though erroneous T 0 , K, and ρc are obtained at the second calculation stage, the goodness of the fit and the extrapolated temperatures can be improved substantially.Examples on synthetic as well as on field data demonstrate the accuracy of the algorithm.With the use of the algorithm, we can determine the ambient equilibrium temperature and extrapolated temperatures to within 0.1 mK from synthetic data and within a few 0.1 mK from field data which is better than that provided by Hartmann and Villinger (2002).Particular emphasis is laid on the effect that the solutions for these variables obtained from five-digit resolution (0.1 mK, high-resolution) data exhibits a great increase in accuracy over that from four-digit resolution (1 mK) data.
These four parameters can be substituted into Eq. 1 for extrapolation beyond the used temperature recording and has been tested in a field experiment.Comparing the extrapolated temperature at 2.38 hours after the penetration with the recorded data, the discrepancy is about 0.1 m°C.In deed, the extrapolation can also be used to remove the residual heat that may affect the accuracy of calculating conductivity during the heat-pulse experiment (heat-pulse decay).
) It may converge slowly.(2) It may oscillate widely, i.e., continually change directions.(3) It may not converge at all.In order to overcome the shortcomings, we propose a two-stage solution process, that is, to solve T B , T 0 and K first with the empirical relationship between K and ρc ,

Fig. 3 .
Fig. 3. Comparison of the extrapolated temperatures ( , × , ) from different resolution levels of data with the given model ( − ).Note the difference is dramatically reduced by using the high-resolution data (0.1 mK).

Fig. 4 .
Fig. 4. Comparison of the extrapolated temperatures ( , , × , +) from different time-shifts data with the model ( − ).It demonstrates that the tempera- tures can be accurately extrapolated (less than 0.1mK in error) by the proposed method and are almost not affected by the time-shifts.

Fig. 5 .
Fig. 5.A photograph of Lister-type heat flow probe used in the experiment.Inside the tube, thermistor spacing is about 60~100 cm (depending on the length of the tube).Total weight is about 350 kg.

Fig. 6 .
Fig. 6.Comparison of the recorded temperature data measured with a typical resolution (1 mK) heat probe (top) with those measured with a highresolution (0.1 mK) heat probe (bottom).For emphasizing the difference, the temperature scale and the deviation with respect to the regression lines have been greatly enlarged.

Fig. 7 .
Fig. 7. Comparison of the extrapolated temperatures with the recorded temperatures.Recorded temperature data ( ) at 1460 m water depth were obtained during a heat flow measurement cruise off southwestern Taiwan.The extrapolated data ( ) regressed using the frictional heat raised temperature data ( ) are very close to those of field data .As expected, the extrapolated temperatures asymptotically approach the solved equilibrium ambient temperatures T B1 and T B 2 ( − ) (therefore the gradients) of the sediment.