Numerical Modeling of the Interference of Thermally Unbalanced Aquifer Thermal Energy Storage Systems in Brussels (Belgium)

Abstract

A numerical model was built using FEFLOW^® to simulate groundwater flow and heat transport in a confined aquifer in Brussels where two Aquifer Thermal Energy Storage (ATES) systems were installed. These systems are operating in adjacent buildings and exploit the same aquifer made up of mixed sandy and silty sublayers. The model was calibrated for groundwater flow and partially for heat transport. Several scenarios were considered to determine if the two ATES systems were interfering. The results showed that a significant imbalance between the injection of warm and cold water in the first installed ATES system led to the occurrence of a heat plume spreading more and more over the years. This plume eventually reached the cold wells of the same installation. The temperature, therefore, increased in warm and cold wells and the efficiency of the building’s cooling system decreased. When the second ATES system began to be operational, the simulated results showed that, even if the heat plumes of the two systems had come into contact, the influence of the second system on the first one was negligible during the first two years of joint operation. For a longer modeled period, simulated results pointed out that the joint operation of the two ATES systems was not adapted to balance, in the long term, the quantity of warm and cold water injected in the aquifer. The groundwater temperature would rise inexorably in the warm and cold wells of both systems. The heat plumes would spread more and more over the years at the expense of the efficiency of both systems, especially concerning building’s cooling with stored cold groundwater.

1. Introduction

New demands for renewable energy sources have greatly increased the attention given to shallow geothermal systems. When groundwater pumping and reinjection in a shallow aquifer are considered associated with the use of a heat pump, an Aquifer Thermal Energy Storage (ATES) can be developed. The basic operation of an ATES system is the following: in summer, groundwater is pumped for cooling. Heat is captured from the building and transferred to water which is then injected back into the aquifer at a higher temperature. In winter, the system works the other way. Groundwater is pumped from the aquifer to heat the building and injected back into the aquifer at a lower temperature. The advantage of ATES stands in the seasonal storage of heat and cold around warm and cold wells, respectively. However, very specific hydrogeological conditions are still required [1,2]: the aquifer hydraulic conductivity should be high enough to allow significant pumping to provide the required power, the groundwater flow should be limited enough to prevent the migration of the heat plume far away from where it can be used. For large power demands (P > 100 kW), the efficiency and costs of an ATES system are known as quite optimized as compared to Borehole Thermal Energy Storage systems known also as ‘closed-loop’ geothermal systems. While these systems are currently favored in urban areas to heat and cool big buildings, leading in turn to CO₂ savings [3,4,5], permit considerations could be a major obstacle in rural areas as the protection of groundwater resources is often mostly influenced by drinking water regulations [1].

The renewability and efficiency of ATES systems could be threatened by inadequate conceptual and operational conditions. The objective of this paper is to show how the thermal characteristics of the ATES system of a building and the hydrogeological conditions can interact. Here, we will see that the remaining heat plumes due to the poor thermal balance of the building (i.e., which requires more cooling than heating yearly), are not removed by a sufficiently important groundwater flow in the aquifer. This is probably instructive to show such a kind of ‘worst case study’ as it induces an unsatisfactory situation impacting further ATES systems (even if these latter are balanced). Thus, the illustrated case study is showing how one thermally unbalanced ATES system could partly jeopardize shallow groundwater conditions in an urban area where other systems could suffer from the induced thermal pollution limiting both efficiency and durability.

2. Materials and Methods

2.1. Context

Nowadays, more than 2800 ATES systems are operating worldwide [3]. In total, 99% of ATES systems are low-temperature ATES and are working at temperatures below 25 °C [3]. Most of them are located in Europe, and particularly in the Netherlands [6]. In Belgium, 20 ATES systems are already used in Flanders in 2011 [7], four are operating in Brussels [8], and the first one is currently under construction in Wallonia [9]. The number of ATES systems is relatively low in Belgium, mainly because of the conflictual use of aquifers for ATES and drinking water supply, raising concerns among drinking water companies and environmental regulators [1,7,10].

Specific legislations are still relatively vague or non-existing in many countries including European countries where the EU water directive (EU-WFD) mentioned that the release of heat in the underground could be considered as groundwater pollution [11,12]. Only a few countries such as Germany, Denmark, Sweden, and Switzerland have national regulations notably about temperature limits and/or minimum distances between systems, varying widely [11,12]. As mentioned by Bloemendal et al. [13], the only ‘coordination’ takes place ‘through government issuing permits on a first-come, first-served basis’.

An elaborated harmonized management framework structure and a governance model, which provide a roadmap for the different levels of management development of 13 EU countries was achieved in 2020 in the framework of the MUSE project [14]. This sustainable governance model is adaptable to the different urban scales and independent of hydrogeological conditions. The generalized structure of the shallow geothermal energy management framework adopted in MUSE allows the effective analysis of policies to identify potential problems and plan effective solutions, as well as to select the best management objectives, strategies, and measures according to the proposed policy principles.

Even if ATES allows better use of energy and a reduction of CO₂ emissions and has proven its cost-effectiveness [3,6], uncertainty about the long-term effects of thermal storage in the subsurface can hinder the development of this technology [15,16]. Indeed, the difference between the temperature of the injected water and the undisturbed temperature of the aquifer leads to the formation of thermally affected zones. These cold and hot plumes can be modeled by using analytical solutions [17,18] and coupled groundwater flow and heat transport numerical models [19,20,21,22,23,24,25,26,27,28,29,30].

Environmental impacts related to the induced cold or warm zones in aquifers are not yet fully understood [31,32] especially in terms of groundwater quality and biodiversity changes affecting various ecosystems [33]. Furthermore, thermal plumes can reach other wells, leading to thermal interferences. These plumes can be observed inside a system or between adjacent systems and can have a negative or a positive impact on the performance of the systems, whether the interfering wells are storing, respectively, water at different or similar temperatures [4,34]. Negative interferences should therefore be avoided choosing proper well-to-well distances and adapted pumping rates [35,36] to the local hydrogeological conditions. Thermal interferences can also result from thermal imbalances, induced by unequal yearly volumes of warm and cold water injected in the aquifer [37]. Cooling demand may exceed heating demand even in Central and Northern European countries [38]. Such variations in the heating and cooling demand are indeed related to climatic fluctuations but are also induced by inefficient management of the heating system [39] or inadequate structural and architectural features of the building.

Some research has been undertaken to find out how to perform adequate planning to maximize the number of systems without creating negative interferences between them [4]. Indeed, a large distance between wells is better to avoid interferences and to increase the performance of the system but reduces the number of systems that can be installed in a given area [4,34]. It seems also that only 50–60% of the authorized ATES capacities are used on average in the Netherlands [13,39]. Oversizing is applied mainly because of the uncertainty related to the predicted energy demand of buildings during their lifetime [13].

As ATES is a recent technology, there is a lack of long-term monitoring data of operating ATES systems [6,39]. This lack is induced by the high costs of measuring devices and the difficult maintenance of permanent monitoring wells in urban areas [5]. Additionally, the optimized use of the thermal storage and heating/cooling system to provide the maximum comfort of the building occupants does not correspond automatically to the optimized use of the subsurface storage hydrogeological conditions [13,39]. Ideally, networks of ATES systems could also be studied and optimized towards an optimum for both the subsurface and comfort conditions in the buildings.

2.2. Study Case

2.2.1. Problem Statement

Two ATES systems exploiting the same confined aquifer and located in adjacent buildings have been operating for a few years. The objective is to evaluate the potential thermal interferences between the two systems, using a groundwater and heat transport numerical model with the available operational data. The first system was installed in Building n° 1 (Figure 1) and started to operate in March 2014, while the second one was installed in Building n° 2 (Figure 1) and started to operate around August 2017. The two systems are characterized by eight wells each. Four of them pump groundwater to heat the building in winter, heat being also stored around these wells seasonally (‘warm well’) by injection of warm water in summer. The four other wells are used to pump groundwater to cool the building in summer, cold being also stored around these wells seasonally (‘cold wells’) by injection of cold water in winter.

The numerical model was implemented in FEFLOW^®, this software allows to simulate groundwater flow in saturated or variably saturated conditions, coupled with the transport of heat and contaminants. Using the finite element numerical method, this code allows a highly flexible meshing strategy and solves the equations described hereunder (see Section 2.3.1) with temperature-related water density and viscosity variations. The available operational data are total pumping rates measured during the operation of the ATES system in Building n° 1 and a few temperature measurements from the wells of Building n° 1, representative of the temperature of the aquifer. For the second building, the power and energy derived from groundwater have been recorded. From these data, groundwater pumping and injection rates have been deduced (i.e., using Equation (4)).

2.2.2. Hydrogeology Context

The study area is located in the Brussels-Capital Region of Belgium, characterized by marine Cenozoic deposits overlying Mesozoic Cretaceous chalk and the Paleozoic bedrock [40]. The shallow hydrogeological context of Brussels consists, therefore, in a succession of sand and clay deposits (

Table 1. Geology as observed during drilling at the location of Building n° 1.

Altitude (m DNG)	Depth	Description	Formation
+18 to +11	0 to 7	Loam, clay	Quaternary deposits
+11 to −2	7 to 20	Sand with gravel	Quaternary deposits
−2 to −5	20 to 23	Clayey sand to sandy clay	Lower Eocene
−5 to −40	23 to 58	Clay	Lower Eocene
−40 to −46	58 to 64	Sand	Upper Palaeocene
−46 to −56	64 to 74	Sand and clay	Upper Palaeocene
−56 to −67.5	74 to 85.5	Clay	Upper Palaeocene
−67.5 to −76.5	85.5 to 94.5	Marl (above) and chalk (below)	Cretaceous
From −76.5	From 94.5	Mainly phyllites with quartzitic phyllites intercalations	Paleozoic basement

), i.e., alternating aquifers and aquitards overlying the chalk aquifer and the Paleozoic fractured aquifer [41]. A simplified map showing the location of the studied site and the boundaries of the model (see Section 2.3.3) is shown in Figure 2 (the local background map of Brussels is not shown due to confidentiality constraints). The two ATES systems are exploiting the Cenozoic Palaeocene aquifer system divided into multiple sub-layers with the main fine sand layers separated by a clayey layer (Table 1). The sandy parts act as aquifers and host the well screens, while the central unit is an irregular clayey aquitard. Groundwater interactions with the Cretaceous chalk and the Paleozoic bedrock are considered negligible based on the very low permeability of the lower layer of the Upper Palaeocene (aquiclude) and observations of quite independent piezometric behaviors in these aquifer systems compared to the upper system [42]. Furthermore, locally, pumping tests were performed in a well exploiting the upper aquifer system without any measured groundwater arrival and influence from the Cretaceous chalk and the Paleozoic bedrock aquifer systems.

2.3. Methods

2.3.1. Groundwater Flow and Heat Transport Equations

The 3D transient groundwater flow in saturated conditions can be described by the following equation [1]:

$$\nabla \cdot \rho \left(\mathit{K}\cdot \nabla h\right)+\rho q^{\prime }=\rho S_s\frac{\partial h}{\partial t}$$

(1)

where $\mathit{K}$ is the hydraulic conductivity tensor (m/s), $\nabla h$ is the piezometric gradient vector (-), $\rho$ is the groundwater density (kg/m³), $q^{\prime }$ is the water flow rate per unit volume of the geological medium (s⁻¹) that is injected ($q^{\prime }$ > 0) or withdrawn ($q^{\prime }$ < 0), $S_s$ is the specific storage coefficient (m⁻¹) which is the volume of groundwater gained or released per unit volume of a saturated porous medium and per unit of h change.

The heat transport equation can be written [1]:

$$\frac{\partial {\rho }_b{c}_{b}T}{\partial t}=-\nabla \cdot \left[{\rho }_w{c}_w\mathit{q}T-{\lambda }_b\nabla T\right]+Q_T$$

(2)

where $T$ is the temperature (°K) and $\nabla T$ the temperature gradient (°K/m), ${\lambda }_b$ is the heat (or thermal) conductivity (W/(m°K)) of the bulk porous medium, ${\rho }_w$ is the water density (kg/m³), $c_w$ is the water heat capacity (J/(kg°K)), so that ${\rho }_w{c}_w$ is the water volumetric heat capacity (J/(m³°K)), $\mathit{q}$ is the total water flux vector (m/s) from Darcy’s law and the possible temperature effect on the water density and viscosity, ${\rho }_b$ is the bulk density (kg/m³) of the porous medium, $c_b$ is the bulk heat capacity (J/(kg°K)) of the porous medium, so that ${\rho }_b{c}_b$ is the bulk volumetric heat capacity (J/(m³°K)), $Q_T$ is the heat source (if $Q_T$> 0) or sink (if $Q_T$< 0) term (W/m³).

Indeed, coupling effects between these flow and heat transport equations must be taken into account through the dependency of the hydraulic conductivity values on the water viscosity (µ) and density variations induced by temperature changes:

$$\mathit{K}\left(T\right)=\mathit{k}\rho \left(T\right)g/\mu \left(T\right)$$

(3)

When producing thermal energy with a heat pump, the power P (in W) can be expressed by [1]:

$$P=\frac{Q∆T{\rho }_w{c}_w}{\left(1-\frac{1}{COP}\right)}$$

(4)

where COP is the coefficient of performance (-) of the heat pump (i.e., the ratio between the useful produced thermal work and the required work), Q is the water flow rate in the heat pump (m3/s), ∆T the temperature difference between the upstream and downstream of the heat pump (°K). In practice, the power P is usually expressed in kW.

2.3.2. Conceptual Model

The groundwater model boundaries were defined regarding the hydrogeological context and considering the boundary conditions (BCs) that could arise from them. The best option was to place boundaries along with physical barriers or the furthest possible from the place where results are expected, to avoid boundary conditions influencing these results [1]. In this study case, the top and the bottom of the model correspond to the top and bottom, respectively, of the Cenozoic aquifer as this is confined between two very low permeability confining clayey layers. For the lateral boundaries, since the aquifer is confined over a large area, an actual change of hydrogeological conditions could not be advocated to define boundaries. These boundaries were therefore set at a ‘large enough’ distance from the study area, which was determined with regard to the measured radius of influence during previous pumping tests performed in the site of interest. The lateral boundaries were defined largely outside the measured drawdown cones induced by these tests. Since the tests were performed without injection and with larger pumping rates than those observed during the ATES operations, we can consider that the lateral boundaries are defined conservatively. The SE and NW lateral boundaries (Figure 2) are prescribed piezometric head BCs set along piezometric head contour lines of the sand aquifer, as delineated in previous projects in Brussels [42]. Through the NE and SW lateral boundaries (Figure 2), zero-flux BCs are considered as they are set along with perpendicular directions with regard to the main regional groundwater flux. Then, at the location of the ATES pumping and injection wells but also the location of the other groundwater abstraction wells located inside the modeled zone, prescribed fluxes are imposed, corresponding to the averaged values of pumping/reinjection during the different considered periods.

Concerning the heat transport simulations, heat BCs must also be prescribed. The available temperature data over the modeled zone being limited, a temperature was prescribed along all lateral boundaries corresponding to the undisturbed background temperature in the aquifer. This latter was measured at the beginning of the operation of the first ATES system and is equal to 12.4 °C. On the upper and lower boundaries, prescribed zero-conduction fluxes are added to zero-advection conditions to result in zero-heat flux BCs. Any influence of the air temperature and from the mean thermal energy earthflow (estimated between 0.05 and 0.11 W/m² [43]) is neglected here.

2.3.3. Numerical Model

The finite element method was used via the FEFLOW software. A finite element mesh was designed, defining nodes at locations of ATES wells and abstraction wells. The mesh was refined around all wells and a higher refinement level was used around the zone of interest, leading to a mesh of about 62,000 elements (Figure 3). The model is 3D extended to a third dimension with 10 layers representing the Cenozoic aquifer system. Based on detailed drilling information from the different boreholes, the first 5 layers are corresponding to the top part of the aquifer, layers 6 and 7 correspond to the clayey part and layers 8–10 to the lower aquifer part. The aquifer system is fully confined and therefore always saturated.

Groundwater flow simulations were first performed to calibrate the model. Steady-state simulations were carried out to obtain the initial hydraulic heads for transient state simulations. Once the model was calibrated for groundwater flow, heat transport was taken into account to simulate the operation of the two ATES systems and to obtain the study results.

The boundary conditions defined for groundwater flow simulations are (Figure 2):

• on the upper and lower boundaries: zero-flux BCs;
• along the SE and NW boundaries: imposed piezometric heads, representative of the piezometric level observed during the simulated period;
• along the NE and SW boundaries: zero-flux BCs;
• in abstraction wells: imposed flux, equal to the average annual abstracted flow rate.
• The defined heat transport BCs are (Figure 2):
• on the upper and lower boundaries: zero-diffusive flux BCs;
• along lateral boundaries: imposed temperature, equal to the undisturbed background temperature.

Steady-state conditions were first calculated from initial confined piezometric conditions. The calculated results were then used as initial conditions for the transient state simulations. In simulations involving heat transport, the undisturbed initial background temperature of 12.4 °C was chosen over the initial temperature in the whole domain.

3. Results

3.1. Groundwater Flow Calibration

The model was first calibrated for groundwater flow. Values of the input parameters were adapted manually to obtain simulated heads as close as possible to piezometric heads observed in monitoring wells in the area of interest (Figure 2). Calibration was first performed in steady-state to obtain initial conditions for the transient state calibration. An automatic calibration was also used in a steady-state but not in a transient state. Two different periods were selected to calibrate the model in a transient state. These periods correspond to pumping test periods having led to important drawdowns in the considered aquifer. To simplify the calibration process in a transient state, only the input parameters (hydraulic conductivity and specific storage coefficient) of the first five layers were modified. In the clayey part (layers 6 and 7) and in the lower aquifer part (layers 8–10) typical values equal to 1 × 10⁻⁸ and 1.6 × 10⁻⁵ m/s, respectively, were adopted for the hydraulic conductivity and 5 × 10⁻⁵/m for the specific storage coefficient [42].

Two examples of results are given in Figure 4 showing for the two piezometers Pz 367 and Pz 518 (Figure 2), located the closest to the considered ATES systems, the transient evolution of the computed versus measured piezometric heads during the two computed periods of 2013 and 2016.

The hydraulic conductivity values obtained in the five upper layers are ranging between 2 × 10⁻⁶ and 3 × 10⁻⁵ m/s and the specific storage coefficient values are ranging between 1 × 10⁻⁵ and 5 × 10⁻⁴/m (Figure 5). With such parameter values, the obtained simulated piezometric heads were the closest to the observed ones.

3.2. Heat Transport Simulations

Once the model was calibrated for groundwater flow, heat transport was added to simulate the operation of ATES systems in Building n° 1 and Building n° 2.

The ‘Open-Loop Design plug-in’ was used in FEFLOW to simulate the operation of ATES systems. Power demands and injection temperatures, therefore, had to be defined. For Building n° 1, the power extracted from the subsurface was calculated based on abstraction rates measured during the operation of the system and assuming a ΔT of 6°K between extraction and injection temperatures. For Building n° 2, the power was determined from the energy derived from the subsurface and assuming a ΔT of 6°K as well. This assumption was made on injection temperatures because the temperatures measured during the operation of the second ATES system were not available and those measured during the operation of the first ATES system were not reliable. Indeed, measurements being performed at the entrance and exit of the heat exchanger, are influenced by the ambient temperature of the building basement and by frequent switches between groundwater injection and abstraction. For Building n° 1, the obtained power curves already highlight an imbalance between cooling and heating powers (Figure 6). The cooling season is globally longer than the heating season and the power demand is way larger for cooling than for heating. This imbalance is not observed in the power demand curve of Building n° 2 (Figure 7), where the power demands are comparable for heating and cooling and the duration of the seasons is also nearly equivalent. These differences are probably due to huge structural and architectural differences between the two buildings.

Due to a lack of temperature observation data, calibrating the model for heat transport was limited to check that the simulated temperature values and trends are realistic when compared to observed values (Figure 8) for the reference scenario (see hereunder). Typical heat transport parameter values were introduced in the model (

Table 2. Heat transport parameters values introduced in the model [1,2,7,8].

Hydrogeological Unit	Lithology	Layers	${\mathit{n}}_{\mathit{e}}$(-) (Effective Transport Porosity)	${\mathit{\rho }}_{\mathit{s}}{\mathit{c}}_{\mathit{s}}$ (MJ/m3/°K) (Solid Heat Volumetric Capacity)	${\mathit{\lambda }}_{\mathit{s}}$(W/m/°K) (Solid Heat Conductivity)
Upper aquifer part	Sand	1–5	0.1	2	2.5
Low K part	Clay	6–7	0.03	3	2
Lower aquifer part	Sand	8–10	0.1	3	2

).

3.2.1. Scenario Descriptions

Different scenarios were considered for the heat transport simulations. The first six scenarios only take the operation of the first ATES system into account. The operation of the second ATES system was added in scenarios 7 and 8. In scenario 1, taken as the reference scenario, the operation of the first ATES system was simulated from 1 January 2014 to 29 February 2020. To investigate the sensitivity of the model results to changes in values of the heat transport parameters, scenarios 2, 3, and 4 were performed. In scenario 2, the thermal conductivity of the solid matrix was increased from 2.5 to 4 W/m/°K in layers 1–5 and from 2 to 2.5 W/m/°K in layers 6–10 of the model. In scenario 3, the volumetric heat capacity of the solid matrix was modified from 2 to 3 MJ/m³/°K in layers 1–5. Additionally, in scenario 4, the effective transport porosity was increased from 0.10 to 0.15 in layers 1–5.

Then, taking back the heat transport parameters from scenario 1 (reference scenario) the effect on the results of the ΔT used to calculate the power was investigated in scenarios 5 and 6, using a ΔT of 4 and 8°K, respectively. In these two scenarios, the power is the same as in the reference scenario, but the flow rates are different since these latter depend on the abstraction and injection temperatures.

Then, taking back the condition from scenario 1 for the first ATES system, operation of the second ATES system, starting in August 2017, was added in scenario 7, to detect whether or not this second system has an impact on the first one and if the two systems interfere. Finally, a longer scenario, scenario 8, was considered to see how the thermal plumes would behave in the future if current operating conditions are maintained. This last simulation was performed until 31 December 2029, with average powers obtained from the previous operation of the systems and a ΔT of 6°K in both cases.

3.2.2. Simulation Results

In the first scenario (reference scenario), the temperature increases year after year in warm (WW) and cold wells (CW) (Figure 8). The same trend could be observed in the measured temperatures. Starting from the natural background temperature of groundwater around 12.5 °C, a step increase was simulated in the hot wells at the spring period (May), corresponding to the start of the building cooling period and thus to injection of hot water in the hot wells. Then, the temperature in the hot wells is increased during the whole summer period. During this summer period, the temperature is also increased in the pumping cold wells to converge back towards the natural background groundwater temperature. This latter is true from 2014 until 2016, but in the following years, the coldest spring temperatures in the cold wells start to be shifted higher and higher, and thus, during the summer, the temperatures increase towards higher temperatures each year. That is clearly due to the Building n° 1 thermal imbalance that requires more cooling than heating, especially during the last years. Similarly, at the starting of the winter period (October), a step temperature decrease is logically simulated in the cold wells corresponding to the start of the building heating period and thus injection of cold water in the cold wells. Then, the temperature in the cold wells is decreased during the whole winter period. During this winter period, the temperature is also decreased in the pumping hot wells to converge back towards the natural background groundwater temperature but never reaching the initial natural temperature due to the thermal imbalance. Consequently, and year after year, the coldest spring temperatures in the hot wells are increased from 12.5 to 18 °C.

During the operation of the ATES system, heat and cold were stored through the injection of warm and cold water. However, the amount of warm water injected into the subsurface was larger than the amount of cold water. This thermal imbalance led to thermal interferences between the heat plume and the cold wells of this system (Figure 9) and the temperature increased in the cold wells. The cold storage is therefore ineffective or non-existent and the efficiency of the system for cooling decreased year after year.

The modifications of the thermal conductivity made in scenario 2 had practically no effect on the temperatures observed in the various wells. The heat plumes observed in the two cases at similar periods are identical (Figure 9). When a higher volumetric heat capacity of the solid matrix was used (scenario 3), a slight decrease in the temperature could be observed compared to the reference scenario. Indeed, as the volumetric heat capacity increases, more energy is needed to heat the solid matrix and less thermal energy reaches the cold wells. The simulated heat plume at the end of the last simulated season is slightly less spread out in scenario 3 than in scenario 1 (Figure 6) but the influence on the results remains however limited. With a higher effective transport porosity in layers 1–5 (scenario 4), the temperatures observed in the wells of Building n° 1 are the same as those in the reference scenario. The changes in heat transport parameters had almost no effect on the simulation results (Figure 6).

With a lower ΔT to calculate the power derived from the subsurface (scenario 5), the amplitude of the temperature variations between the heating and cooling seasons is lower. The temperature increases more slowly in the warm wells, but the heat plume still reaches the cold wells (Figure 9). Inversely, with a higher ΔT (Scenario 6), the amplitude of the temperature variations is higher than in the reference scenario. Higher temperatures are observed in the heat plume reaching the cold wells (Figure 9). However, the temperatures obtained with a ΔT of 6°K are closer to measured temperatures than the temperatures observed in the other two cases.

When the operation of the second ATES system is taken into account (scenario 7), the temperatures observed in the wells of Building n° 1 stay similar to those in the reference scenario. During the first two years of its operation, the second ATES system had therefore almost no effect on the operation of the first ATES system. However, both systems interact as their heat plumes came into contact in 2019. The second system is more balanced than the first one since heat and cold plumes are induced (Figure 9). However, in the long term, the results of scenario 8 show that, with the current operational conditions, the temperature will increase in the warm and cold wells of both systems. The heat plumes will enlarge year after year (Figure 9), at the expense of cold storage and efficiency for cooling of both systems. Even if the second system is more balanced, the cooling season is still 15 days longer than the heating season, which could lead to a larger injection of warm water than cold water.

Since the flow is from southeast to northwest, the hot water plume created by Building n° 1 (in the northwest zone) could be dissipated further in the northwest direction. Unfortunately, the groundwater flow is relatively weak, and the hot plume does not migrate significantly. On the other hand, the hot water plume created by Building n° 2 could potentially affect partially the cold wells of Building 1. In all scenarios, the simulated plumes keep a relatively circular shape showing the limited influence of the groundwater flow. Thermal energy losses by advection are therefore limited. In this particular case, heat is not transported away from the wells, and the temperature increases locally around these wells. As the groundwater flow influence is quite limited, it would be even more necessary to reach a thermal balance between abstraction and injection. If it is not the case, heat or cold accumulates and finally reaches the opposite wells. The system, therefore, becomes inefficient.

4. Discussion

The heat transport simulations that were performed indicated that the large imbalance between the injection of warm and cold water in the subsurface led to thermal interferences between the heat plume and the cold wells of Building n° 1. This resulted in a temperature increase in the warm and cold wells of this system. The temperature in cold wells is now higher than the undisturbed background temperature in the aquifer, which leads to a lower efficiency of the system. Changing the main heat transport parameter values did not influence significantly the results, which demonstrated that even if these values were slightly different in reality, this would not have a large impact on the results. The simulations made with a different ΔT showed that this parameter has a large impact on the results and that the temperatures obtained with a ΔT of 6°K were the closest from the measured temperatures. Adding the operation of the second ATES system had almost no effect on the temperatures and heat plumes of Building n° 1. The two systems, however, interacted as their heat plumes came into contact. In the long term, keeping operational conditions the same as currently, the heat plumes will enlarge year after year, at the expense of cold storage and efficiency for cooling of both systems. The local hydrogeological conditions are such that the growing heat plume is not counterbalanced by sufficient groundwater flow to bring the groundwater temperature down to its background value.

We are clearly showing here a bad example, but we believe that it is by pointing out such examples that we can avoid this type of mistake in the future. Here, the management of both ATES systems should therefore be changed quickly to rebalance the heat and cold storage. This change must be designed to limit interferences between wells of the same system and/or of the neighbor system and to improve efficiency in the short and long terms. This is especially true for Building n° 1 since the ATES system in Building n° 2 has been in operation for less than three years and no imbalance in heat and cold storage has been observed yet.

Indeed, for Building n° 1, this is particularly uneasy to correct such a thermal unbalance after construction and active ATES operations for years. To balance the amount of heat stored in the aquifer, less warm water and/or more cold water should be injected. Thus, in other words, the use of the system for cooling should be decreased, and/or the use of the system for heating should be maximized. Among the possible solutions currently investigated for Building n° 1, there are:

• increased use of heat during winter by providing heat to neighboring buildings;
• ncreased night ventilation and window solar protection during the summer.

These efforts should reduce the need for cooling during warm days combined with an increased need for heat during winters. At the same time, new aquifer temperature measurements (temperature sensors) are being installed to obtain better monitoring of the impact of the ATES systems. Different management scenarios involving the above-mentioned solutions should then be simulated for prediction purposes and further decisions.

More generally, if the concept of seasonal heat storage is considered in an aquifer, results from this bad case study are striking by the induced impact in terms of heat plume. Here, the groundwater flow and heat advection are very limited, and a global thermal balance for hot and cold injections would have been required before issuing the permitting certificate. On the contrary, in other hydrogeological conditions, if the groundwater flow and heat advection are very important, thermal energy cannot be stored efficiently locally as heat and cold plumes are transported far away from the wells [36,44,45,46]. If the groundwater flow and heat advection are moderate, a detailed simulation of the groundwater flow and heat transport in the aquifer is particularly required to find out if the annual imbalance can be managed in relation to the specific local hydrogeological conditions.

Another important point lies in the fact that the reliability and accuracy of the model results are indeed highly related to the available measurements: piezometric heads, groundwater temperatures, detailed pumping, and injection flow rate, etc. After a few months or years, these data can be used to improve the reliability of the model results through a longer and more detailed calibration procedure. This last point is often forgotten or ignored by real estate developers and civil engineering designers in their rush to build as quickly as possible and at the lowest cost.