Estimation of the vertical borehole thermal parameters based on the evolution algorithm using temperature response functions

The vertical borehole ground heat exchange performance is still the issue for the engineers installing vertical borehole ground plants using ground-coupled heat pumps. Besides geological and climate change challenges, they face the extraordinary heat transfer process between the circulating fluid and the ground surrounding the U-tube and interactions of vertical boreholes. This paper describes the technique to evaluate the grout, soil thermal parameters and borehole thermal resistance simultaneously using the particle swarm optimization algorithm. The reference thermal response data set from the sandbox laboratory was used for the analysis. A thermal response test (TRT) was made, including the different temperature response functions, for a few time scales. The estimates and errors of the grout, soil thermal conductivity and borehole resistivity were presented and compared with the results of the laboratory experiment and researchers. The target functions, in our case root mean square error values, were less than 0.034 for all analysis cases. The calculation algorithm was written using the Matlab 2016 program and could be easily expanded by increasing the number of target functions and evaluation algorithms. The presented TRT data analysis will increase the knowledge about the vertical borehole ground heat exchange design.


INTRODUCTION
It is well known that numerical methods are more accurate than analytical methods for the design and simulation of the heat transfer of borehole ground heat exchangers because the heat transfer rates varying in time and the influence of the surrounding ground are taken into account.Despite the lower accuracy and higher computational costs, the analytical methods are still popular for in situ thermal response test analysis.The combinations of analytical and numerical methods are usually implemented into borehole ground heat exchanger (BGHE) design and optimization software.The most popular programs are Earth Energy Designer (EED) (Eskilson 1987), GLHEPRO (Spitler 2017), TRNSYS (Klein & Beckman 2007), EnergyPlus (Crawley et al. 2001).The temperature response function, well known as a 'G-function', was first mentioned by Eskilson (1987).The G-functions could be separated into time periods for the short (Yazuzturk & Spitler 1999;Zeng et al. 2002;Lamarche & Beauchamp 2007), intermediate (Carslaw & Jaeger 1986) and large (Ingersoll 1954;Eskilson 1987)  , and   ,   and  are ground thermal diffusivity, borehole radius and active length, respectively.The finite line source models (Zeng et al. 2002;Claesson & Javed 2011;Li et al. 2014) are most applicable and efficient in modelling the underground heat transfer process.Li et al. (2014) and Claesson & Javed (2011) developed a heat transfer model covering time scales from minutes to decades.Bandos et al. (2009) described the Gfunctions for the mid-point of a borehole  =  2 and the average temperature along the BGHE z-axis for intermediate and long time scales.Practical applications are directly related with the publications by Popov et al. (2016) and Zhang et al. (2017).
The results of the present study were partly presented by the authors at the annual Lithuanian Mathematical Society conference in 2015.The research presents a new technique for the calculation of the thermal parameters using the particle swarm optimization algorithm and testing the validity of borehole heat transfer analytical algorithms.The introduction gives an overview of previous implementation results.The methodology part presents different time scale G-functions and the stochastic evolution algorithm.The experiment section provides the TRT data analysis, including the geological environment having the reference data set from sandbox laboratory experiments.The TRT reference data set was used for the simultaneous calculation of the soil, grout thermal conductivity and borehole thermal resistivity using the evolution algorithm.The estimates and errors of parameters are provided in the section of results, which were compared with laboratory thermal probe results and data from other publications.The conclusion and remarks are presented in the last section, together with the further steps for investigation and research.

Heat transfer model using temperature response functions
The theoretical heat transfer model was proposed by Li et al. (2014) for  unit-step heat load where   is the mean temperature of circulating inlet and outlet fluid in the U-tube,  0 is undisturbed field temperature,   is the borehole thermal resistance,   () is the temperature response function for the infinite line source model.The borehole thermal resistance was calculated by the formula (Hellström 1991) where  = (2 •   /  ) 4 ,   ,   are the grout and the U-tube pipe thermal conductivities, respectively; ℎ  is the convective heat transfer coefficient of fluid,   ,   and   denote the borehole inner radius, and the outer and the inner radius of the U-shaped pipe,   is the spacing between the centre of the legs of the U-type tube.The heat transfer models incorporating different temperature response G-functions will be described below.

Infinite line source (ILS) model
Starting from Carslaw & Jaeger (1986) an analytical approach was developed for mean temperature on the borehole wall.This equation is most commonly used for the estimation of the thermal conductivity of the ground during in situ TRTs.The temperature response function   () for  unit-step heat load is derived by the formula where   is the thermal diffusivity of the surrounding ground,   is the thermal conductivity of the surrounding ground,  is time and  is the integral variable.The   (t) suffers from the impact of ground surface temperature variation in the thermal process for a long time period.

Infinite cylinder source (ICS) model
The infinite cylinder with a fixed heat flux rate should be used for BGHE approximation where the U-pipe could be approximated as 'equivalent diameter'.The ILS method mentioned above is a simplified version of the ICS method.The temperature response function can be described following Ingersoll (1954): where J0, Y0, J1, Y1 are Bessel functions of the first and second kind, and  = are the G-function parameters.More analytical solutions were defined by Carslaw & Jaeger (1986); other authors (Kavanaugh & Rafferty 1997) discussed the design of BGHEs.Claesson & Javed (2011) presented the following formulation for the integral mean temperature at a distance r of a finite length line heat source (FLS) extending from  =  to  =  + .The surface  = 0 is maintained at a temperature  = 0.The G-function at the distance  =   of the borehole at time t has the expression

Finite line source (FLS) model
where () denotes a complementary function and   () is a temperature response function as the average temperature on the borehole wall.The transient thermal process between the ground surrounding the borehole and the backfilling material in the borehole could not be accounted for in equations ( 4), ( 5) and ( 6).
The G-functions estimates for giving the maximum errors of the time criterion   ≥ 5  2   do not exceed 10%.Gehlin (2002) showed that the maximum error could be less than 2.5% for   ≥ 20  2   .The estimates of the average fluid temperature are calculated as shown in formula (8): where G(t) can have one of the formulations mentioned above:   () (4),   () (5),   () (6).

Evolution algorithm: particle swarm optimization
In the year 1995 the stochastic optimization technique particle swarm optimization (PSO) was proposed by Kennedy & Eberhart (1995).The introduced algorithms imitated the social behaviour of a flock of birds.
Typically, the population of particles is called a swarm which consists of  particles moving in a problem search space.Each particle is defined as a potential solution.For an N-dimensional search space, the position of the th particle is represented as    = ( 1 ,  2, … ,   ).At each generation, the new particle position is found by adding a displacement to the current position where the displacement could be calculated by equation ( 9) where    and   +1 represent the current and previous positions of particle ,   +1 is the current velocity of particle  and is represented as    = ( 1 ,  2 , … ,   ).For each generation the velocity of each particle is updated by the formula where    and   +1 are the current and previous velocities of each particle , and inertial weight was changed for every iteration as  =   − Here  ,  is the average fluid temperature during the experiment,  ,  is the estimate of fluid temperature enabling heat transfer simulation,  is the number of TRT test data points.Equation (10) consists of three parts: the first one is called the momentum part that defines the previous velocity, the second part is called the cognition part that represents the best position of an individual particle, the third part is called the social component that represents the collaboration among particles in the swarm.The cognitive learning coefficient is   and the social learning coefficient is   , and  1 ,  2 are two random numbers generated by uniform distribution within the interval [0, 1].The relative sizes of these components determine their contribution to the new particle velocity.It is well known that the standard PSO algorithm can be balanced between the global and the local minimum because of the proper selection of the inertial weight parameter.Clerc (1999) suggested how to assure the convergence of the algorithm for the determination of heat transfer coefficients.The aim of our analysis is to iteratively estimate the unknown heat transfer coefficients using the PSO procedure, which results in a negligible difference between temperature measurements taken at the given locations and temperatures computed from the numerical simulation.The numerical simulation temperatures were calculated on the borehole wall using the abovementioned thermal response functions.The fitness function RMSE value of each particle at the th iteration is given by the difference between the measured and calculated temperature curves, at the position    .A short description of the PSO algorithm is given in Fig. 1.

EXPERIMENT: IN SITU TRT SETUP
Many thermal response tests are performed in real in situ geological conditions.The high quality of reference data can assure the quality of results.The reference data set could be used for testing and validating the new transfer model.We used the reference data set from a large laboratory sandbox for testing and validating the heat transfer parameters.The sandbox was constructed from a wooden frame.The sandbox is a rectangle with sides of 1.83 m and 18.32 m.The borehole was centred horizontally along the length (18.32 m) of the sandbox.A plastic liner separated the sand from the wooden frame in order to keep water.The sand was saturated from the local utility water line by five perforated parallel lines uniformly spaced on the bottom of the wooden box.The external parts of the wooden sandbox were wholly thermally insulated to minimize the effect of changing weather conditions.The parameters of the borehole installed into the sandbox are shown in Table 1.
The TRT test was designed so that the heat input rate and the fluid flow rate through the closed U-pipe were close to a constant value.The U-pipe was installed into the horizontal borehole (tube) made of aluminium and the distance between the U-pipe centres was fixed following the high quality of TRT test procedures.Table 1 shows the technical parameters of the experimental TRT apparatus that were used for the TRT test.The reference data set consists of 24 thermistors that provide temperature measurements at the borehole wall and at specific locations in the surrounding soil, and the inlet (21) and outlet (20) of fluid (Fig. 2).Measurements are recorded every minute on a computer data acquisition system for the heat transfer model.Beier et al. (2011) described the measurement procedures determining grout and ground thermal conductivities by using a non-steady-state thermal probe invented by Hooper & Lepper (1950).The estimated uncertainty was ±5% for grout and ground thermal conductivities.The same TRT reference data set was used by Javed (2012) for the analysis and validation of the borehole heat transfer model.
A testing unit for in situ thermal response tests is connected to the U-tube in the sandbox.Two electric heating elements together supply approximately 1056 W to the circulating fluid during the TRT test.The pipe material, grout and thermal properties are given in Table 2.
The bentonite grout having 20% solids was mixed with water in order to make the borehole filling.Before each TRT test, the uniform temperature of air, fluid and ground should be measured.The fluid circulation through the U-tube should start together with the heating by electric elements, which ensure a constant heat input rate to the water.All measurements of the temperature at the locations of thermistors, fluid flow rate and heat input were connected and recorded in a data acquisition system once per minute.The TRT test was conducted for about 52 h within the saturated sandbox and serves as a reference data set for heat transfer modelling and simulation.

RESULTS
We applied the standard PSO algorithm to the TRT data set of the single-objective minimization problem in the two-dimensional search space composed of the ground and grout thermal conductivity parameters.The number of optimization problems was equal to the different temperature response functions   (),   (),   () which were used in the heat transfer simulation.All these benchmark functions together with the objective RMSE function were performed in the Matlab 2016 programming language.The numerical simulation was performed using these parameters as follows: the number of particles N = 100; inertia weight  ∈ (0.2; 1.2); the particle and swarm best weights   ,   = 2, respectively.The maximum number of generations and simulation runs was equal to 20.The mean values of the surrounding air temperature and heat flux per unit of borehole length were used for the performance of heat transfer simulation.Before data analysis, the grout and ground thermal conductivity values in the two-dimensional search space could vary in the interval [0;5].Some practical investigations were done to define the inertia weighting function, and particle and swarm weight values in advance.First, the efficiency of the linear decreasing inertia weight formula shown by Bansal et al. (2011) was used in the PSO algorithm.Second, the trial and error method which gives good results but is not always the rule of thumb was used to select the particle and swarm weights.We carried out the linear independence analysis of thermal parameters, the impact of the U-pipe shank spacing value on the effective borehole thermal resistance and calculated errors of thermal parameters under various TRT durations.Before the estimation of the thermal parameters the linear dependence analysis should be performed and analysed.The first derivatives of fluid temperature with respect to thermal parameters were defined in order to get the relative sensitivity coefficients (RSCs) defined by formula ( 12) The Matlab 2016a Symbolic toolbox was used to calculate the RSCs.In these calculations the assumptions of linear independence analysis should be satisfied.First, there is no possibility of writing RSCs of thermal parameters making a linear combination of other RSCs.Second, relatively small RSC values can lead to |  | ≈ 0 and the thermal parameters cannot be estimated simultaneously (Zhang et  2018).The relatively large value of the determinant can assure the linear independence in estimating the thermal parameters (Ozisik 2018).The RSC values for heat transfer models using different G-functions   ,   and   were calculated and presented in Fig. 3A-C.The high |  | values (see Fig. 4) assure the linear independence between the thermal parameters and the PSO algorithm can be performed by selecting them simultaneously.The exact RSC values of grout thermal conductivity and U-pipe shank spacing are -7.261 and -7.735, respectively.
The durations of TRT reference datasets were divided by selecting different starting points 1h, 2h, 7h, etc. excluding the first 2 h from analysis.Gehlin (2002) made many practical investigations concerning TRT test duration.In practice the ILS model is most appreciated but the TRT experiment should not be less than 50 h and the first few hours are excluded from analysis.
The TRT experiment duration is equal to  ≤ 5  2   , which was proposed by Zhang et al. (2018).The authors validated this proposal by showing the importance of the TRT data starting-point and the stability of thermal parameter estimates.The errors are less than 3% for the mean value of the identified soil thermal parameters if the TRT duration is not less than 28 h.
The estimates of thermal parameters using heat transfer model G-functions for different TRT startingpoints and durations are illustrated in Fig. 5A-C.
The uncertainty analysis was performed to get necessary knowledge about the borehole thermal resistance due to the uncertainty of the installed U-pipe location.The U-pipe shank spacing was identified using Hellström (1991) formula.The designed U-pipe spacing between centres   (0.053 m) was changed to 0.0688 m.The estimated and calculated   values are shown in Fig. 6.The proposed value of   is suitable for performing the heat transfer modelling and eliminates the uncertainty of the fluid temperature prediction.
The statistics of thermal parameter estimates for a stable period of more than 28 h are given in Table 3.The errors are benchmarked with Beier et al. (2011) and Zhang et al. (2018) results.The mean RMSE values are around the uncertainty of temperature measurement.
In the investigations by Beier & Smith (2003) and Beier et al. (2011) relative errors around 5% are considered as reliable.Zhang et al. (2018) presented the applicability of genetic algorithms despite the high relative errors.

DISCUSSION
The soil and grout parameters were estimated by the PSO algorithm using TRT reference set data.The new proposed method was performed for three heat transfer models   (),   (),   () using various TRT test durations with different starting-points.
The new method for the evaluation of thermal parameters provides acceptable relative errors for soil, grout thermal conductivities and borehole thermal resistivity simultaneously.The TRT duration should not be less than 28 h, which was directly related with the TRT data starting-point.The PSO algorithm did not fall into the local minimum and showed a good performance in identifying the thermal properties of parameters of porous media.The research conducted has certain real significance, using the actual TRT data from different time intervals.The relative error of borehole thermal resistance is under 4%, which is affected by U-pipe shank spacing and soil and grout thermal parameters.
time scales.The time thresholds conform with the duration of the thermal response test (TRT), which can vary from 40 to 240 h.The time shifts are  <   for short time,   ≤  ≤   for intermediate time and  >   for long time.Here

Fig. 1 .
Fig. 1.Estimation of the thermal parameters by the PSO procedure.

Fig. 3 .
Fig. 3. (A) Relative sensitivity coefficients of thermal parameters using the GILS function.(B) Relative sensitivity coefficients of thermal parameters using the GICS function.(C) Relative sensitivity coefficients of thermal parameters using the GFLS function.

Fig. 4 .
Fig. 4. Determinant of relative sensitivity coefficients under various heat transfer methods.

Fig. 5 .Fig. 6 .
Fig. 5. (A) Estimates of thermal parameters using the GILS method for different TRT durations.(B) Estimates of thermal parameters using the GICS method for different TRT durations.(C) Estimates of thermal parameters using the GFLS method for different TRT durations.
, with   and   being inertial maximum and minimum values, respectively.The previous best position of each particle could be defined as  , , giving the best fitness function value.The global best position   = ( ,1 ,  ,2 , … ,  , ) is

Table 1 .
Description of the sandbox experiment

Table 2 .
Experimental values of thermal parameters 2 /s Fig. 2. Schema of measurement points on one half of the sandbox. al.

Table 3 .
Statistics of identified thermal parameters