Performance assessment of a simplified hybrid model for a vertical ground heat exchanger

Energy and Buildings 66 (2013) 437–444

Contents lists available at ScienceDirect

Energy and Buildings journal homepage: www.elsevier.com/locate/enbuild

Performance assessment of a simplified hybrid model for a vertical ground heat exchanger Ismael R. Maestre, F. Javier González Gallero ∗ , Pascual Álvarez Gómez, J. Daniel Mena Baladés Escuela Politécnica Superior de Algeciras, University of Cádiz, Avenida Ramón Puyol, s/n, Algeciras 11202, Spain

a r t i c l e

i n f o

Article history: Received 1 April 2013 Accepted 12 July 2013 Keywords: Ground heat exchanger Thermal resistance and capacity model 1-D hybrid model g-Function Building simulation program

a b s t r a c t A simplified model to simulate single U-tube ground heat exchangers is presented in this paper. The model is based on the electrical analogy to simulate heat transfer within the borehole and the use of thermal response factors (g-functions) to estimate injection–extraction heat flow to surrounding ground. The inclusion of two capacities for the heat carrier fluid and the grout material makes the model suitable for short time step simulations. The thermal resistance of the borehole has been changed considering different heat conduction paths. Several variants of the model have been validated through comparison with a refined computational fluid dynamic (CFD) reference model. Good results have been obtained for the different variables analysed (outlet temperature and surface borehole temperature). RMSE (Root Mean Squared Error) values were smaller than 1 ◦ C while relative errors were under 5%. © 2013 Elsevier B.V. All rights reserved.

1. Introduction In 2008, the European Union put in place the Climate and Energy Package to combat climate change, whose targets are known as ‘20-20-20’: 20% reduction in greenhouse gas emissions, raising the share of EU energy consumption produced from renewable sources to 20% and 20% improvement in the EU’s energy efficiency (commission staff working paper [1]). In this context, the ground source heat pump (GSHP) systems, using borehole heat exchangers (BHEs), are renewable energy systems promoted by the Renewable Energy Directive (Directive 2009/28/EC of 23 April 2009, 2009 OJ L 140/16 [2]) that can offer higher energy efficiency for air conditioning compared to conventional air conditioning (A/C) systems. GSHP systems have three main components: the ground side where the BHE is placed and that is used to get heat out of or into the ground, the heat pump to convert that heat to a suitable temperature level, and the building side transferring the heat or cold into the rooms. A heat pump is used in winter to extract heat from the relatively warmer ground and pump it into the conditioned space. The process may be reversed in summer, extracting the heat from the conditioned space and sending it out to a ground heat exchanger that warms the relatively cool ground (Florides and Kalogirou [3]). Due to the high thermal inertia of the soil, ground temperature below a certain depth remains nearly constant throughout the year. Thus, the earth offers a nearly steady and large heat source, heat

∗ Corresponding author. Tel.: +34 956 028000; fax: +34 956 028001. E-mail address: [email protected] (F.J. González Gallero). 0378-7788/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.enbuild.2013.07.041

sink and heat storage medium for thermal energy uses while air temperature is affected by higher climatic variations. More importantly, the use of the soil as a heat source or sink decreases the difference between the reservoir temperatures, getting higher efficiencies. Efficiency of a GSHP system is measured by the coefficient of performance (COP), which can be defined as the ratio between the useful heat output and the electric power input used to drive the compressor and pump for the ground loop. The heating COPs typically range from 3.5 to 5 (Omer [4]) while greater cooling COPs have also been reported (Aikins and Choi [5]). Furthermore, ground heat exchanger can substitute the refrigeration towers avoiding problems related to legionellosis. Blum et al. [6] also reported that the use of GSHP systems in comparison to conventional heating systems leaded to CO2 savings between 15% and 77% depending on the supplied energy for the heat pumps and the efficiency of installation. Thus, GSHP have received an increasing interest in North America and also in the EU, being one of top regions in the development of this kind of technology, with more than 720,000 units at the end of 2007 and an installed power of 8758 Mw (EurObserv’ER Report, [7]). A simultaneous modelling of the whole system (ground side, heat pump and building site) is required to estimate the COP of the overall system. Thus, for a fast and reliable simulation of the system, it is essential to have sufficiently accurate and computationally cheap simulation models of the BHE. The BHE usually consists of a closed loop or coaxial pipes lowered into a borehole. Pipes in a ‘U’ loop and grouted into vertical boreholes are likely the most common form of BHEs. Their design

438

I.R. Maestre et al. / Energy and Buildings 66 (2013) 437–444

Nomenclature c C d f Fo g H h k l L ˙ m Nu Pr Q R r Re t ts T V x

specific heat (J/(kg K)) thermal capacitance (J/K) spacing between pipes (m) friction factor Fourier number adimensional response factor (g-function) borehole depth (m) convective heat transfer coefficient (W/(m2 K)) thermal conductivity (W/(m K)) effective conduction path pipe length (m) fluid flow rate (kg/s) Nusselt number Prandtl number heat flux (W) thermal resistance (K/W) radius (m) Reynolds number time (s) time scale temperature (K) volume (m3 ) adimensional variable equal to ln(t/ts )

Greek symbols ˛ thermal diffusivity (m2 /s) ˇ adimensional factor t discretization time step (s) limiting time scale (s)   density (kg/m3 ) Subscripts b borehole fluid f g grout G ground in inflow index pipe p P pressure Superscripts i time step index time step index k

is critical to the long-term performance of the heat pump system and the application of dynamic models is required to capture the heat transfer inside and outside the borehole. This heat transfer process must be treated on the whole as a transient process. Generally, the main objective of the BHE thermal analysis is to determine the temperature of the heat carrier circulating fluid under different operation conditions. Several models have been developed to simulate ground heat exchangers: analytical, numerical and hybrid models. A detailed literature review of the research, developments and the most typical simulation models of the vertical GCHP and BHEs is presented by Florides and Kalogirou [3] and Yang et al. [8]. Analytical models are based on a number of simplifying assumptions of the borehole and are applied to both the design of BHEs and the analysis of in situ test data (He et al. [9]). The cylinder source model (Carslaw and Jaeger [10]) and further simplifications, such as the line source model (Ingersoll et al. [11]) use heat transfer theory of these sources (or sinks) to solve for the heat transfer rate from

the borehole wall to the surrounding soil neglecting the internal region of the borehole. These simplified models are not suitable for taking account of time varying heat transfer rates and the influence of surrounding boreholes on long time scales. Numerical models and totally discretized models (Bauer et al. [12], Al-Khoury et al. I [13], Al-Khoury et al. II [14]) use finite difference or finite-volume methods to solve for temperature distribution in the whole domain. These kinds of models get accurate simulation results also on short time scales, but they usually require a great computational effort and are not suitable to be implemented in energy building simulation Software program. Many of the models currently used are hybrid ones. Interesting approaches were presented by Eskilson [15] and later by Yavuzturk and Spitler [16], which used dimensionless temperature response factors, the so-called g-functions, for the numerical estimation of the time-dependant borehole outer wall temperature. Thermal resistance and capacity (RC) models (Bauer et al. [12], Huber and Schuler [17]), based on the electrical analogy, can also be an interesting alternative because they combine the advantages of the former models: simplicity and capability to analyze long and short time behaviour. Within this kind of models, those which lump the pipes together into one single centred pipe (of equivalent diameter) oversimplify the complex heat transfer process that takes place inside the grout material, especially because the effect of important parameters like the shank spacing (or distance between pipes) is not considered (Sharqawy et al. [18]). However, although this effect has been analyzed for steady-state conditions, no results have been reported yet for the error made by the use of models with one thermal capacity node (in the fluid) in operating conditions of relatively fast changes in inlet fluid temperature. The present work tries to quantify the error made by the use of a simple BHE model, based on the RC discretization of the borehole, and the use of short-term and long-term g-functions to model ground heat conduction. Model validation is carried out by comparison of its output with the results obtained by a three-dimensional CFD numerical model.

2. Model description This section describes the one-dimensional thermal (radial) model developed for the simulation of BHEs with single-U-pipes. The model uses the RC formulation for the heat transfer within the borehole, with a capacity node both in the fluid and inside the borehole. Then, their thermal inertia is taken into account, which allows the calculation of the BHE response to fast changes of the heat injection–extraction rates. Thermal response of the ground, which is assumed to be homogeneous, is assessed by using short time step and long time-step g-function models. The main hypotheses of the model presented here are described next. The non-homogeneous soil surrounding the borehole is replaced by homogeneous soil with mean transport properties such as mean conductivity and mean diffusivity. These mean values are derived from experimental data. To make the problem onedimensional, the two pipes are replaced by a single pipe with an effective radius. Then, vertical heat conduction can be neglected within the whole simulation domain and the symmetry of the model configuration allows to solve the one-dimensional energy equation in cylindrical coordinates. Boundary condition in the ground is an imposed temperature boundary condition (TG ), equal to the average vertical ground temperature. The superposition principle has been used to consider variable heat injection–extraction rates. The model developed in this work is a simplified and computationally cheap model that, as it will be assessed later on, can be

I.R. Maestre et al. / Energy and Buildings 66 (2013) 437–444

439

Fig. 2. Description of the RC model with one capacity node in the fluid and one in the grout material of the GHE.

where h is the thermal heat transfer coefficient and Af is the (convective) thermal transfer surface of the real pipe. The heat flux in the pipe is conductive (qcd p ) and the corresponding thermal resistance (Rp ) is estimated as the thermal resistance of a cylinder with a length L in steady-state radial heat conduction (Eq. (2)). That is:

Fig. 1. Real configuration of the BHE (with single U-pipes).

considered suitable for the simulation of the short-time step and long-time step response of BHEs (with single U-pipes). Furthermore, the model can be considered versatile because any kind of discrete input is allowed (in terms of both temperature and mass flow), requiring no previous signal processing. 2.1. BHE characterization Fig. 1 shows longitudinal and transverse sections of the BHE, together with some of its geometric characteristics. The real configuration of the BHE is mainly defined by the following geometric parameters: total length of the pipe (L), which is equal to 2H, approximately, where H is the borehole depth; inner and outer pipe radius (rp1 and rp2 , respectively); distance between the pipes (d) and the borehole radius (rb ). The RC model proposed for the BHE and its main physical parameters are shown in Fig. 2. Five thermal capacity nodes in total have been considered in order to characterize the mean fluid temperature (T¯ f ), inner and outer pipe surface temperatures (Tp1 y Tp2 ), grout temperature (Tg ), outer surface temperature of the borehole (Tb ) and the far-field ground temperature (TG ). Mean fluid temperature (T¯ f ) has been defined as the average between the input and the output fluid temperature. Furthermore, Fig. 2 shows the different thermal resistances in each zone of the heat exchanger (fluid, pipe and grout material) and the equivalent thermal capacities of the fluid (Cf ) and grout material (Cg ). The parameters that define the proposed simplified model are ˙ and fluid volume (Vf ) are described next in detail. Mass flow (m) considered equal to those in the real configuration of the BHE. Heat flow (qcf v ) and the thermal resistance (Rf ) within the core of the fluid, which are mainly convective, can be estimated as follows (Eq. (1)): qcf v =

T¯ f − Tp1 Rf

,

Rf =

1 hAf

(1)

qcd p =

Tp1 − Tp2 , Rp

Rp =

ln(rp2 /rp1 ) 2kp L

(2)

Finally, thermal resistances within the grout material (Rg1 and Rg2 ) depend on the effective length of the heat conduction path from the pipe to the borehole surface (lb ), where node b is located (Tb ). The value of the radius rg , which defines the position of the temperature node of the grout material (Tg ), is estimated considering the equality of the volumes of grout material at both sides of the node. Then, for a fixed value of lb , the position of node g will be given by the following equation:

 rg =

2 lb2 + rp2

(3)

2

The value of lb will vary between two limit values that depend on the effective heat conduction path within the inner zone of the BHE. These extreme heat conduction paths have been displayed in Fig. 3. From the corresponding value of rg to each lb value (Eq. (3)), heat flux (qcd ) and thermal resistances Rg1 and Rg2 can easily be g1 calculated (Eq. (4)): qcd g1 =

Tp2 − Tg , Rg1

Rg1 =

ln(rg /rp2 ) , 2kg L

Rg2 =

ln(lb /rg ) 2kg L

(4)

For the estimation of their extreme values, the corresponding limit values of lb and rg must be used in Eq. (4). The equivalent grout capacity is calculated as follows (Eq. (5)): 2 Cg = ˇ · Cg real = ˇg cpg (rb2 − 2rp2 )H

(5)

where ˇ is a user-defined factor that allows the modification of thermal inertia. Cgreal represents the real thermal capacity of the grout material.

440

I.R. Maestre et al. / Energy and Buildings 66 (2013) 437–444

Fig. 3. lb values for different heat conduction paths within the borehole.

2.2. Energy balance in the fluid

Heat flux from node g to borehole surface (node b) is expressed by Eq. (11):

Considering the possibility of energy accumulation in the fluid and a finite difference discretization scheme, the equation of energy conservation leads to the following expression (Eq. (6)): cPf V

T¯ fi − T¯ fi−1 t

=−

(6)

=

+

t



T¯ i − Tgi



 f





 

1/ 2rp1 Lh + ln rp2 /rp1 / 2kp L + ln rg /rp2 / 2kg L



i i ˙ Pf Tin ˙ Pf 2T¯ fi − Tin +mc − mc





(7)

The heat transfer process within the fluid is mainly convective. The convective heat transfer coefficient (h) can be evaluated through correlation equations that take the Reynolds number into account (De Carli et al. [19]). In the model proposed here, the convective heat transfer coefficient is calculated by Eq. (4) for turbulent flows (Eq. (8)): Nu =

g cPg Vg

T¯ fi − T¯ fi−1



(f/8) · Re · Pr 1.07 + 12.7 · (f/8)

1/2

· (Pr 2/3 − 1)

(104 < Re < 5 × 106 )

(8)

The friction factor f is the Darcy friction factor, and Petukhov’s formula (Petukhov [20] Eq. (5)) can be used for the estimation: f = [0.790 · ln(Re) − 1.640]−2

(11)

Rg2

Thus, Eq. (12) can be obtained within the grout material:

i i ˙ Pf Tin ˙ Pf (2T¯ fi − Tin = −qcf v + mc − mc )

The superscript i defines the time when a variable is estimated. Eq. (6) can be expressed in terms of the temperature at the thermal capacity node g, considering the equality of the thermal fluxes defined by Eqs. (1), (2) and (4). Thus, Eq. (7) is obtained: cPf V

Tgi − Tbi

cd qcd g2 = qb =

2.3. Energy balance in the borehole The second thermal capacity node of the model is located within the grout material, allowing the consideration of its thermal inertia, which is especially useful during fast changes of the boundary conditions.

t







T¯ i − Tgi

 f





 

1/ 2rp1 Lh + ln rp2 /rp1 / 2kp L + ln rg /rp2 / 2kg L



Tbi − Tgi

 

ln lb /rg / 2kg L





(12)

The second term of Eq. (12) describes the mean heat flux at the borehole surface (qk ), for which the temperature at the outer wall of the borehole (Tb ) is evaluated at a former time Tbi−1 , instead of at the current time (Tbi ).

2.4. Ground modelling Two kinds of g-functions have been used in the model proposed here: one for small time-steps and the other for long time-steps. Long time-step response factor model developed by Eskilson [15] has been used to estimate g-functions for longer periods. The two asymptotic approximations proposed for the g-function are the following:

glt

t r  b

(9)

During the stop intervals in which the flow speed is nil, h is evaluated using the free-convection correlation for vertical cylinders.

Tgi − Tgi−1

, ts H

⎧  H  1 t ⎪ ⎨ ln 2r + 2 ln ts , b =   ⎪ H ⎩ ln , t > ts

5rb2 /˛ < t < ts (13)

2rb

Consequently, this g-function model is valid for changes of the heat injection–extraction rates that take place after the time limit defined by Eq. (14): =

5rb2 ˛g

(14)

For a borehole in typical applications this time step might be of several hours (De Carli et al. [19]).

I.R. Maestre et al. / Energy and Buildings 66 (2013) 437–444

In order to assess the response of the system to faster changes (t < ), the g-function used in the present work is described by Eq. (15): gst

t r  b , ts H

= ax3 + bx2 + cx + d

(15)

where x = ln(t/ts ) and a, b and c are constants that must be determined. The numerical approach described by Yavuzturk and Spitler [16] can be used for the estimation of the short time-step gfunction. 2.5. Resolution process Using the superposition principle, the BHE temperature response to any arbitrary heat injection–rejection function to the ground can be calculated. In the model presented here, the summation has been splited into two terms, considering if the time interval is higher or smaller than the limiting time scale , and using the corresponding approach of long time-step (glt ) or short time-step (gst ) g-function, respectively. This procedure is described in Eq. (16): Tbi

= TG +

i  qk − qk−1

2kG

k=i−i

 qk − qk−1 k=1

2kG

gst



i−i −1

+



glt

t i − t k−1 rb , ts H

t i − t k−1 rb , ts H

441

Table 1 Configuration and thermal properties of the Ground heat exchanger (GHE). Properties

Value

rb , borehole radius (mm) rG , ground radius (mm) rp1 , pipe inner radius (mm) rp2 , pipe outer radius (mm) d, spacing between pipes (mm) Fluid kf , conductivity (W/m K) cP , thermal capacity (MJ/m3 K) Pipe kp , conductivity (W/m K) p cPp , thermal capacity (MJ/m3 K) Grout kg , conductivity (W/m K) g cPg , thermal capacity (MJ/m3 K) Ground kG , conductivity (W/m K) G cPG , thermal capacity (MJ/m3 K) ˙ (kg/s) Fluid flow rate, m Convective coefficient, h (W/m2 K) Initial ground temperature, TG (◦ C)

75.0 1500–18,000 13.1 16.0 32.2 0.61 4.17 0.39 1.77 0.75 3.90 2.50 2.50 0.54 4949.6 10.0



(16)

The subscript i defines the time at which the time interval equals , that is, t i − t i = . Initial grout temperature (Tg0 ) and borehole temperature (Tb0 ) are assumed to be equal to the ground temperature (TG ). With these initial conditions, Eqs. (7) and (12), which approximately describe the heat transfer processes within the borehole, can be easily solved for each time “i”, leading to the values of T¯ fi and Tgi . Afterwards, Tbi value is estimated by using (16), continuing the marching process. 3. Validation methodology This section describes the validation process of the onedimensional model developed in the present work. This process has mainly consisted of the comparison of the simplified model results with those obtained by a reference model developed using CFD. This reference model, which has a high level of detail and a high computational cost, was done using ANSYS CFX 13.0 software package. The response of the proposed model for short term and long term periods and at the most frequent operating conditions of this kind of systems was analyzed. The output variables under analysis have been the fluid output temperature and the borehole surface temperature. 3.1. Study cases The case under study is a typical domestic building BHE (He et al. [9]). It is a single borehole with a diameter of 150 mm and a depth of 100 m. Geometric and thermal properties of the borehole have been shown in Table 1. The spacing between (d) pipes is the end-to-end distance (see Fig. 1). Although in some study cases the time variation of the mass flow has been considered, its nominal value has been set to 0.54 kg/s, leading to an average fluid speed of 1 m/s. Under these conditions, Reynolds number (Re) is about 26,000, leading to a turbulent flow regime. As the inlet fluid temperature is equal to 20 ◦ C, the system has been assumed to operate in cooling mode.

For the configuration chosen, the limiting time scale  has a value of 7.8 h, approximately. The validation tests proposed here take into account the variation of both the mass flow and the inlet fluid temperature, allowing the simulation of the usual operating conditions of this kind of devices, with stops during the night. A complete validation process, consisted of three different tests, i ) has been carried out. In test No. 1, the inlet fluid temperature (Tin has been kept constant. In test No. 2, a variable inlet mass flow has been considered, using a square step function that ranges from zero to a fixed mass flow of 0.54 kg/s, and a period of 12 h. In test No. 3, a rectangular pulse function has been used, with the same range of mass flow and time intervals of 8 h (zero flux) and 16 h (non-zero flux of 0.54 kg/s). A set of simulations of the BHE has been run with total simulation times of 48 h and 1 year, for all the tests formerly described. The discretization time interval was set to 5 min. 3.2. Reference model The reference model used for the validation of the simplified model has been the CFD model of the BHE configuration described in Section 3.1 of the present work. The overall extension of the soil was chosen big enough to avoid the influence of the outermost element boundary condition. The surrounding ground has been considered of cylindrical shape with a radius of 1.5 m (rG ), for simulations of 48 h and 18 m for 1 year simulations. These rG values were obtained imposing Fourier number (Fo) values smaller than 0.1, in order to assure that the ‘far-field’ ground temperature (TG ) would not be disturbed (zero heat flux) during the total simulation period (48 h and 1 year). This condition is defined by Eq. (17), where ˛G is the ground thermal diffusivity. Fo =

˛G · t rG2

≤ 0.1

(17)

The numerical mesh consists of hexahedron and wedge-shaped volume elements, about 135,177 cells in total (see Fig. 4). The refining of the mesh led to changes lower than 0.5%. The minimum size of the elements was set to 1.5 mm to capture the geometry of the pipes properly. The size of the soil elements was chosen fine enough for them not to influence the transient simulation results considerably. 400 layers of 0.25 m have been used to cover the 100 m depth of the system. The mesh has been checked for cell skewness and aspect ratio according to normal CFD standard. During the

442

I.R. Maestre et al. / Energy and Buildings 66 (2013) 437–444

Fig. 5. Short time-step g-function and polynomial fitting.

Fig. 4. Details of the mesh used to represent the borehole.

simulations all residuals were run down to 10−5 . The k-Epsilon model has been used to simulate turbulence within the pipes. Outer surface of the ground cylinder surrounding the BHE was assumed to be at constant temperature of 10 ◦ C. For test No. 1, a constant inlet fluid temperature of 20 ◦ C was set. For tests No. 2 and No. 3 an inlet boundary conditions of variable mass flow were considered according to the pulse functions formerly described.

3.3. Short time-step g-function A key point of the model presented here is that a priori set of short time-step g-functions (gst ) is needed. This issue can be considered a small drawback of the simplified model, because a certain set of gst functions must be available for different BHE configurations (shank spacing, borehole radius, conductivity, etc.). There are different alternatives to construct the gst functions. They can be approximated by using the line source or the cylinder source, to calculate a fluid temperature profile versus time which could be input into Eq. (18) (Young [21]): gst

t l  b , ts H

=

2kG (T¯ f − Rb q − TG ) q

(18)

Another method, which offers higher accuracy than the analytical ones is the finite volume method (Patankar [22]). There are several programs, like GEMS2D (Rees et al. [23]) or the one developed by Yavuzturk and Spitler [16], which have implemented that method successfully. In the present validation process, the method described by Yavuzturk and Spitler [16] has been used. Thus, the borehole temperature response to a basic step pulse has been estimated in the former reference model. Afterwards, the temperature response is converted to a series of non-dimensional temperature response

factors. Then, the expression obtained for the gst function is illustrated by the following equation: y = gst

t r  b , ts H

= −0.0069x3 − 0.2026x2 − 1.4927x − 0.0109 (19)

where x is equal to ln(t/ts ).The gst function and its polynomial fitting have been shown in Fig. 5. 4. Results and discussion Next, the results obtained for the different configurations of the simplified model under the formerly operating conditions, are described and discussed. The error analysis has especially been based on the outlet fluid temperature (Tout ) and the borehole surface temperature (Tb ). Two error measures have been used: root-mean-square error (RMSE) and relative error. The configurations used correspond to the conduction paths described in Section 2.1, considering simulation time periods of 48 h and 1 year. Thus, borehole thermal resistance can take different values: a maximum value and a minimum value corresponding to the longest conduction path (lbmax ) and to the shortest path (lbmin ), respectively, and a mean value equal to the average of the former limit values. Tables 2 and 3, show the RMSE and relative errors made by the use of the simplified model taking the response of the CFD model as a reference. As it can be seen, the variation of the thermal resistance leads to more meaningful changes of the error values than those related to variations of the borehole thermal capacity. Furthermore, it can be inferred that, in all cases, the smallest errors are obtained for the minimum values of the borehole thermal resistance and the maximum value of the borehole thermal capacity. Thus, for that optimum configuration (with minimum thermal borehole resistance and maximum borehole thermal capacity) RMSE values are significantly smaller than 1 ◦ C, for both outlet fluid temperature and borehole surface temperature. Relative errors are also smaller than 5%.

Table 2 RMSE values and relative errors for test No. 1. Simulation time: 48 h. lb

Minimum

Mean

Maximum

Tout

Tb

ˇ = Cg /Cg real

RMSE (◦ C)

εr

RMSE (◦ C)

εr

1.0 0.5 0.0 1.0 0.5 0.0 1.0 0.5 0.0

0.48 0.56 0.69 0.62 0.71 0.83 1.11 0.80 0.90

0.77 0.87 0.97 1.88 1.95 2.03 2.47 2.53 2.60

0.24 0.22 0.32 0.59 0.56 0.59 0.77 0.74 0.76

1.95 1.47 1.23 5.69 5.21 4.81 7.67 7.21 6.79

I.R. Maestre et al. / Energy and Buildings 66 (2013) 437–444

443

Table 3 RMSE values and relative errors for test No. 2. Simulation time: 48 h. lb

Minimum

Mean

Maximum

Tout

Tb

ˇ = Cg /Cg real

RMSE (◦ C)

εr

RMSE (◦ C)

εr

1.0 0.5 0.0 1.0 0.5 0.0 1.0 0.5 0.0

0.48 0.59 0.72 0.67 0.72 0.84 0.76 0.8 0.9

2.37 3.14 4.01 3.58 4.3 5.09 6.69 4.9 5.63

0.23 0.22 0.33 0.49 0.48 0.52 0.63 0.62 0.64

1.63 1.73 2.27 4.37 4.38 4.56 5.8 5.77 5.87

Fig. 6. Isothermal contourlines in a cross section of the borehole at time t = 120 min. Fig. 8. Outlet temperature for test No. 2.

It is worth to note that, despite the simplicity of the model, very small RMSE values are obtained (the highest value is 1.01) for both outlet fluid temperature and borehole temperature, under all the boundary conditions considered in the validation process. The good results achieved by the simplified model when a minimum borehole thermal resistance is considered, can be justified from the analysis of the heat flux through a transverse section of the borehole (see Fig. 6). Thus, during a few hours of simulation time, CFD model reveals that heat flux seems to concentrate along the line that joins the centre of the pipes (which coincides with the shortest path between the pipe and the borehole wall). However, heat flux is less important along the direction normal to the former one. The results obtained by the simplified model for the outlet fluid temperature under the different test conditions have been illustrated in Figs. 7–9. It can be emphasized that a good fitting is achieved even under conditions of fast variations of the inlet temperature. The sensitivity of the model against variations of the borehole thermal capacity has been illustrated in Fig. 10, considering test No. 2 and a minimum borehole resistance. Thus, no

significant variations of the borehole temperature are obtained when borehole capacity is modified. Similar results are achieved for the rest of cases analyzed. Finally, the long-term response of the simplified model is studied, considering an operating configuration identical to that of test No. 3, but with a longer simulation time (1 year) to check the thermal saturation of the ground. In order to assess the goodness of the model fitting, RMSE values have been estimated for different days of the simulation period (days 10, 14, 182 and 365) through comparison with the results obtained by the CFD reference model. The results shown in Table 5 are those of the simplified model for minimum thermal resistance and maximum thermal capacity. As it can be seen, errors are higher than those estimated during the first 48 h (see Table 4), but the response of the model is very suitable during the whole simulation period (RMSE values are smaller than 0.57 ◦ C).

Fig. 7. Outlet temperature for test No. 1.

Fig. 9. Outlet temperature for test No. 3.

444

I.R. Maestre et al. / Energy and Buildings 66 (2013) 437–444

Table 4 RMSE values for test No. 3. Simulation time: 48 h. lb

Minimum

Mean

Maximum

Tout

Tb

ˇ = Cg /Cg real

RMSE (◦ C)

εr

RMSE (◦ C)

εr

1.0 0.5 0.0 1.0 0.5 0.0 1.0 0.5 0.0

0.21 0.23 0.27 0.39 0.41 0.43 0.5 0.51 0.53

2.84 3.78 4.8 4.11 4.85 5.75 4.74 5.39 6.2

0.26 0.19 0.17 0.75 0.69 0.65 1.01 0.96 0.92

1.51 1.73 2.44 3.6 3.75 4.04 4.7 4.8 4.99

References

Fig. 10. Borehole temperature for test No. 2. Influence of the variation of the borehole capacity. Table 5 Assessment of long-range response. RMSE values for test No. 3. Simulation time: 1 year. Day

10

14

182

365

RMSE-Tout (◦ C) RMSE-Tb (◦ C)

0.51 0.34

0.56 0.44

0.40 0.36

0.41 0.40

5. Conclusions A simplified hybrid model has been used to model the thermal behaviour of a BHE with single U-tube. The model is based on a RC strategy for heat transfer within the borehole and the use of gfunctions for the simulation of the extraction-injection heat flux at the ground. Two thermal capacities have been considered to model the response of the fluid and the grout to changes of the boundary conditions. The error of the model has been analyzed through comparison with a high resolution CFD model. Results show good agreement with RMSE and relative error values for out temperatures smaller than 1 ◦ C and 5%, respectively. The simplicity of the model makes it useful for its use in building thermal simulation programs. Acknowledgements This work has been partially supported by the Andalusian Regional Goverment and the European Union through the European Project FEDER-FSE 2007-2013. The authors would like to thank Mrs. Ma Mar González for her help with translation.

[1] Commission Staff Working Paper, Analysis of Options Beyond 20% GHG Emission Reductions: Member State Results, European Commission, Brussels, 1.2.2012. SWD 2012. [2] Directive 2009/28/EC of 23 April 2009, 2009 OJ L 140/16. [3] G. Florides, S. Kalogirou, Ground heat exchangers – a review of systems, models and applications, Renewable Energy 32 (2007) 2461–2478. [4] A.M. Omer, Ground-source heat pumps systems and applications, Renewable and Sustainable Energy Reviews 12 (2008) 344–371. [5] K.A. Aikins, J.M. Choi, Current status of the performance of GSHP (ground source heat pump) units in the Republic of Korea, Energy 47 (2012) 77–82. [6] P. Blum, G. Campillo, W. Münch, T. Kölbel, CO2 savings of ground source heat pump systems – a regional analysis, Renewable Energy 35 (2010) 122–127. [7] The State of Renewable Energies in Europe, 8th EurObserv’ER Report. 2008, EurObserv’ER Barometer, Paris, December 2008. [8] H. Yang, P. Cui, Z. Fang, Vertical-borehole ground-coupled heat pumps: a review of models and systems, Applied Energy 87 (2010) 16–27. [9] M. He, S. Rees, L. Shao, Simulation of a domestic ground source heat pump system using a transient numerical borehole heat exchanger model, in: Eleventh International IBSA Conference, Glasgow, Scotland, July 27–30, 2009. [10] H.S. Carslaw, J.C. Jaeger, Conduction of Heat in Solids, Claremore Press, Oxford, 1947. [11] L.R. Ingersoll, O.J. Zobel, A.C. Ingersoll, Heat Conduction with Engineering, Geological, and Other Applications, McGraw-Hill, New York, 1954. [12] D. Bauer, W. Heidemann, H. Müller-Steinhagen, H.J.G. Diersch, Thermal resistance and capacity models for borehole heat exchangers, International Journal of Energy Research 35 (2011) 312–320. [13] R. Al-Khoury, P.G. Bonnier, R.B.J. Brinkgreve, Efficient finite element formulation for geothermal heating systems. Part I: steady state, International Journal for Numerical Methods in Engineering 63 (2005) 988–1013. [14] R. Al-Khoury, P.G. Bonnier, Efficient finite element formulation for geothermal heating systems. Part II: transient, International Journal for Numerical Methods in Engineering 67 (2006) 725–745. [15] P. Eskilson, Thermal analysis of heat extraction boreholes, in: Department of Mathematical Physics and Building Technology, University of Lund, Sweden, 1987 (Doctoral Thesis). [16] C. Yavuzturk, J.D. Spitler, A short time step response factor model for vertical ground loop heat exchangers, ASHRAE Transactions 105 (2) (1999). [17] A. Huber, O. Schuler, Berechnungsmodul für Erdwärmesonden. Forschungsprogramm Umgebungs-und Abwärme, Wärmekraftkopplung, Bundesamt für Energie, Bern, 1997. [18] M.H. Sharqawy, E.M. Mokheimer, H.M. Badr, Effective pipe-to-borehole thermal resistance for vertical ground heat exchangers, Geothermics 38 (2009) 271–277. [19] M. De Carli, M. Tonon, A. Zarrella, R. Zecchin, A computational capacity resistance model (CaRM) for vertical ground-coupled heat exchangers, Renewable Energy 35 (2010) 1537–1550. [20] B.S. Petukhov, Advances in Heat Transfer, vol. 6, Academic Press, New York, 1970. [21] T.R. Young, Development, Verification, and Design Analysis of the Borehole Fluid Thermal Mass Model for Approximating Short Term Borehole Thermal Response, Oklahoma State University, 2001 (Doctoral Thesis). [22] S.V. Patankar, Numerical Heat Transfer and Fluid Flow, Hemisphere, New York, NY, 1980. [23] S.J. Rees, J.D. Spitler, X. Xiao, Transient analysis of snow-melting system performance, ASHRAE Transactions 108 (2) (2002) 406–423.