61. Interacting Borehole Heat Exchangers in a Geothermal Setting

Shallow groundwater geothermal investigations often include more than one borehole heat exchanger (BHE) (Al-Khoury et al., 2021). In such applications, understanding the thermal interaction among multiple BHEs as well as on the flowing groundwater is made easier with a numerical groundwater flow and heat transport model. In this example, the accuracy of the groundwater energy transport (GWE) model is demonstrated for a convective-conductive porous domain with multiple thermally-interacting BHEs using an analytical solution first published in (Al-Khoury et al., 2021).

61.1. Example description

For this example nine BHEs are arranged in a 3 \(\times\) 3 configuration with a spacing of 5 \(m\) from each other. The grid extent is 90 \(m \times\) 60 \(m\). Each BHE represents a cyclindrical source of heat with energy being added at a rate of 100 \(\tfrac{W}{m}\) using the energy source loading (ESL) package. In order to better simulate the outward propagation of heat from each BHE, the discretization by vertices (DISV) grid type for both the groundwater flow (GWF) and GWE models was employed (fig Figure 61.1). Grid refinement was added around each BHE (fig Figure 61.2). Grid discretization is coarsened toward the perimeter of the model grid.

../_images/ex-gwe-geotherm-grid.png

Figure 61.1 Configuration of the DISV model grid used to demonstrate the use of GWE in a geothermal transport problem. The original numerical model grid was published in (Al-Khoury et al., 2021). The red box shows the location of Figure 61.2 which provides zoomed-in detail of the model grid in the vicinity of the BHEs.

../_images/ex-gwe-geotherm-grid-inset.png

Figure 61.2 A zoomed-in view showing the refined discretization in close proximity to the BHEs. Cell dimensions are on the scale of centimeters around the perimeter of the BHEs. The diameter of each BHE is 10.0 \(cm\).

In order to test the MODFLOW 6 solution against the published analytical solution, the heat transport model simulates a porous media domain with a porosity of 0.20. Within the GWF model, two constant head (CHD) packages were setup on the left and right sides of the model, respectively, to drive groundwater flow from left to right with a velocity of \(1 \times 10^{-5} \tfrac{m}{s}\) (fig Figure 61.3). The initial temperature throughout the model domain is 0.0 \(^{\circ}C\). Energy is added to the grid cell in the middle of the model domain at a rate of 100 \(\tfrac{W}{m}\). Parameters used for the MODFLOW 6 simulation of the geothermal heat transport problem are shown in Table 61.1.

Table 61.1 Model parameters for example ex-gwe-geotherm.

Parameter

Value

Number of periods in flow model (\(-\))

1

Number of layers (\(-\))

1

Simulation width (\(m\))

60

Simulation length (\(m\))

90

Horizontal hydraulic conductivity (\(m/d\))

1.0

Top of the model (\(m\))

1.0

Bottom of the model (\(m\))

0.0

Porosity (\(-\))

0.2

Length of simulation (\(days\))

50

Initial Temperature (\(^{\circ}C\))

0.0

Advection solution scheme (\(-\))

TVD

Thermal conductivity of water (\(\frac{W}{m \cdot ^{\circ}C}\))

0.56

Thermal conductivity of aquifer material (\(\frac{W}{m \cdot ^{\circ}C}\))

2.50

Density of water (\(kg/m^3\))

1000

Heat capacity of water (\(\frac{J}{kg \cdot ^{\circ}C}\))

4180.0

Density of dry solid aquifer material (\(kg/m^3\))

2650.0

Heat capacity of dry solid aquifer material (\(\frac{J}{kg \cdot ^{\circ}C}\))

900.0

Latent heat of vaporization (\(\frac{J}{kg \cdot ^{\circ}C}\))

2500.0

No mechanical dispersion (\(m^2/day\))

0.0

No transverse dispersivity (\(m^2/day\))

0.0

Starting head (\(m\))

1.00

../_images/ex-gwe-geotherm-head.png

Figure 61.3 The established head gradient drives groundwater flow from left to right with a velocity of \(1 \times 10^{-5} \tfrac{m}{s}\).

61.2. Example Results

Results from the geothermal model run are compared to a published analytical solution 50 days after the start of the simulation (Al-Khoury et al., 2021). Isotemperature contours at 1, 2, 3, 4, 6, and 8 \(^{\circ}C\) provide a visual summary of the match between GWE and the analytical solution (Figure 61.4). Isotemperature contours match particularly well at the lower temperatures (\(\leq 2 ^{\circ}C\)). At temperatures \(>2^{\circ}C\), the simulated temperatures have not advanced as far in the downgradient direction as the analytical solution would suggest.

../_images/ex-gwe-geotherm-temp50days.png

Figure 61.4 Simulated and analytical isotemperature contours for the geothermal example problem. The solid lines correspond to the GWE solution while the dashed line represents the analytical solution published by (Al-Khoury et al., 2021).

61.3. References Cited

Al-Khoury, R., BniLam, N., Arzanfudi, M. M., & Saeid, S. (2021). Analytical model for arbitrarily configured neighboring shallow geothermal installations in the presence of groundwater flow. Geothermics, 93, 102063. https://doi.org/10.1016/j.geothermics.2021.102063

61.4. Jupyter Notebook

The Jupyter notebook used to create the MODFLOW 6 input files for this example and post-process the results is: