Inﬂuence of Groundwater on the Very Shallow Geothermal Potential (vSGP) in the Area of a Large-Scale Geothermal Collector System (LSC)

: The water balance in the very shallow subsurface can be inﬂuenced by capillary rise due to a high groundwater table. Since moisture content is an important factor for the thermal conductivity of soils, this can also have an inﬂuence on the very shallow geothermal potential (vSGP). For this reason, the effect of spatial and seasonal variations in groundwater tables on moisture content in essential depth layers was investigated at a large-scale geothermal collector system (LSC) in Bad Nauheim, Germany. Quasi-one-dimensional simulations using the FEFLOW®ﬁnite-element simulation system were employed to determine site-dependent and seasonally varying moisture contents, from which thermal conductivities were derived. The model setup was previously validated based on recorded moisture contents. The simulations resulted in groundwater-related maximum seasonal and spatial differences in thermal conductivity of 0.14 W/(m · K) in the LSC area. Larger differences of up to 0.21 W/(m · K) resulted for different soil textures at the same depth due to different thermal properties. The results indicate that an efﬁcient design of LSCs requires a sufﬁciently detailed subsurface exploration to account for small-scale variations in grain size distribution and groundwater level.


The Significance of the Very Shallow Geothermal Potential (vSGP) in the Context of the Large-Scale Geothermal Collector System (LSC) in Bad Nauheim
A large-scale geothermal collector system (LSC) was installed in two layers (1.5 m and 2.9 m below ground level (m b.g.l.)) in Bad Nauheim, Hesse, in central Germany [1,2]. The approximately 11,000 m 2 area is used for agriculture or as a horse paddock. The LSC, in combination with a low-temperature district heating and cooling network (5GDHC), supplies a new residential area with cooling and heating energy. Approximately 2.3 GWh/a of source heat is required to heat 400 residential units: 50% of this energy is extracted from the subsurface by the LSC, and 50% is obtained by the 5GDHC system due to low temperature levels and uninsulated piping [1]. The very shallow geothermal potential (vSGP) in the form of thermal conductivity is therefore crucial for the heat extraction of the overall system. The project is scientifically accompanied by the research project "KNW-Opt" (Grant No. 03EN3020).
The potential of very shallow geothermal applications can be described by the systemspecific heat extraction [3] or by the thermal conductivity of soils (λ) [4]. The parameter λ is determined by soil physical parameters such as moisture content (θ), bulk density and grain size distribution [4][5][6][7]. Therefore, the thermal conductivity is independent of the geothermal system used and dependent on the available subsurface. The individual The valley deposits overlay the 2.1 to 4.3 m thick bench gravel aquifer [25]. Due to the varying ground elevations, the distance between the groundwater and surface on site can range from 2.0 to 7.0 m. Since October 2020, groundwater levels have been measured with fluctuations of up to 1.9 m at six monitoring wells around the LSC (Figure 1a). Figure 1b shows the different groundwater depths due to the existing slope, despite the low hydraulic gradient seen in Figure 1a.
Geosciences 2023, 13, x FOR PEER REVIEW 3 of 22 agriculturally influenced silt loam. Three open-end tests [23] in the silty clay loam layer at 1.5 m b.g.l. yielded hydraulic conductivities between 2.0•10 -3 and 2.2•10 -4 m/d. Several infiltration tests with double rings [24] yielded higher values, but these are considered less reliable. The valley deposits overlay the 2.1 to 4.3 m thick bench gravel aquifer [25]. Due to the varying ground elevations, the distance between the groundwater and surface on site can range from 2.0 to 7.0 m. Since October 2020, groundwater levels have been measured with fluctuations of up to 1.9 m at six monitoring wells around the LSC (Figure 1a). Figure  1b shows the different groundwater depths due to the existing slope, despite the low hydraulic gradient seen in Figure 1a.
(a) (b) Figure 1. (a) Test site and location of monitoring wells B1 to B6, measuring points with moisture sensors and exemplary groundwater elevation contour lines from groundwater table measurements taken on April 07, 2021, given in meters above sea level (m a.s.l.), adapted from Rammler et al. [22]; (b) exaggerated schematic cross-section through the collector area with stratigraphic units to illustrate the different distances of the groundwater level (GWL) to the collector layers (LSC) at 1.5 m and 2.9 m b.g.l. The location of the cross-section is shown in Figure 1a. Here, valley deposits include topsoil and agriculturally influenced soil.
Sieve analyses based on DIN EN ISO 17892-4 [26] showed contents of gravel or coarser components of between 20% and 80% for the bench gravel layer. The grain size distributions are too non-uniform to derive hydraulic conductivity [27][28][29][30]. Evaluations of the recovery phase of two short pumping tests gave hydraulic conductivity values of 32 and 58 m/d for the bench gravel.
The boundary with the valley deposits is characterized by a gravelly loam or loamy gravel layer, which is stratigraphically largely assigned to the bench gravel. Below the aquifer, decomposed and locally water-conducting basalt is present [25].

Monitoring of Moisture content, Groundwater Table and Weather Data
The KNW-Opt research project involves intensive monitoring of the LSC [2]. Monitoring data on soil, groundwater and weather conditions were used for this study. (a) Test site and location of monitoring wells B1 to B6, measuring points with moisture sensors and exemplary groundwater elevation contour lines from groundwater table measurements taken on 7 April 2021, given in meters above sea level (m a.s.l.), adapted from Rammler et al. [22]; (b) exaggerated schematic cross-section through the collector area with stratigraphic units to illustrate the different distances of the groundwater level (GWL) to the collector layers (LSC) at 1.5 m and 2.9 m b.g.l. The location of the cross-section is shown in Figure 1a. Here, valley deposits include topsoil and agriculturally influenced soil.
Sieve analyses based on DIN EN ISO 17892-4 [26] showed contents of gravel or coarser components of between 20% and 80% for the bench gravel layer. The grain size distributions are too non-uniform to derive hydraulic conductivity [27][28][29][30]. Evaluations of the recovery phase of two short pumping tests gave hydraulic conductivity values of 32 and 58 m/d for the bench gravel.
The boundary with the valley deposits is characterized by a gravelly loam or loamy gravel layer, which is stratigraphically largely assigned to the bench gravel. Below the aquifer, decomposed and locally water-conducting basalt is present [25].

Monitoring of Moisture Content, Groundwater Table and Weather Data
The KNW-Opt research project involves intensive monitoring of the LSC [2]. Monitoring data on soil, groundwater and weather conditions were used for this study.
Moisture and temperature sensors were installed in the middle of the LSC at 1.5 m, 2.2 m and 2.9 m b.g.l. (measuring point "mLSC") and at the "Reference" measuring point at 0.1 m, 0.25 m, 0.75 m, 1.5 m, 2.2 m and 2.9 m b.g.l. (Figure 1a). For the installation of the moisture sensors, it was necessary to ensure that the sensors were surrounded by fine-grained material. Therefore, to avoid damage, they were embedded in sand.
The area around "Reference" is covered with conifers, resulting in significantly increased water discharge in the shallow subsurface.
Furthermore, pressure and temperature sensors were installed in six groundwatermonitoring wells in the vicinity of the LSC (Figure 1a) with different groundwater distances due to the existing slope of the terrain (Figure 1b). These monitoring wells are described in detail in Bertermann and Rammler [25]. Recorded data are sent from an EASYBus data logger to a central server at 5-min intervals. All measurement components were calibrated by the manufacturer and adjusted to the different cable lengths. The measurement configurations of the used sensors are shown in Table 1. Table 1. Description and data specifications of the sensors used for the moisture content and temperature measurements in the soil and the groundwater level (GWL) and groundwater temperature measurements in the monitoring wells.

Domain
Additionally, groundwater tables were manually recorded almost monthly starting in October 2020 ( Figure 2).
Moisture and temperature sensors were installed in the middle of the LSC at 1.5 m, 2.2 m and 2.9 m b.g.l. (measuring point "mLSC") and at the "Reference" measuring point at 0.1 m, 0.25 m, 0.75 m, 1.5 m, 2.2 m and 2.9 m b.g.l. (Figure 1a). For the installation of the moisture sensors, it was necessary to ensure that the sensors were surrounded by finegrained material. Therefore, to avoid damage, they were embedded in sand.
The area around "Reference" is covered with conifers, resulting in significantly increased water discharge in the shallow subsurface.
Furthermore, pressure and temperature sensors were installed in six groundwatermonitoring wells in the vicinity of the LSC (Figure 1a) with different groundwater distances due to the existing slope of the terrain (Figure 1b). These monitoring wells are described in detail in Bertermann and Rammler [25]. Recorded data are sent from an EASYBus data logger to a central server at 5-minute intervals. All measurement components were calibrated by the manufacturer and adjusted to the different cable lengths. The measurement configurations of the used sensors are shown in Table 1. ±0.5% 0 to 1 bar Additionally, groundwater tables were manually recorded almost monthly starting in October 2020 ( Figure 2). Continuously recorded and manually measured groundwater levels in meters below ground level (m b.g.l.) and their location relative to the LSC layers. According to current knowledge, it is unclear what the temporary short-term stronger fluctuations in B5 can be attributed to. Surface influences could be a possible cause. Also shown are time periods 1 and 2, which were used for the simulations. Period  To monitor weather conditions, a weather station was installed on the roof of the energy center at a height of about 3 m directly east of the LSC. Here, long-wave sky irradiation, short-wave radiation (global radiation), wind direction and speed, ambient temperature, precipitation and relative humidity were recorded at 1-min intervals starting in January 2021. Data gaps and time periods before January 2021 were filled with data from two other climate measurement systems of the German weather service (DWD) delivering open-source data within a range of 15 km from Bad Nauheim (https://opendata.dwd.de/, accessed 27 January 2023).

Quasi-One-Dimensional Simulations with FEFLOW ®
The simulations to investigate the influence of different groundwater levels on moisture content were performed using the FEFLOW ® (v8.0) finite-element simulation system from DHI Wasy GmbH [31]. A quasi-one-dimensional approach was adopted for a detailed investigation at a specific site focusing on vertical water movement. For this purpose, the models were created as one-meter-wide, two-dimensional vertical planar projections with a quadrangle element width of 1 m and height of 0.01 m. Water flow into or out of the model is only possible through the surface or through the gravel aquifer. The schematic model setup is shown in detail in Figure 3.
Geosciences 2023, 13, x FOR PEER REVIEW 5 To monitor weather conditions, a weather station was installed on the roof o energy center at a height of about 3 m directly east of the LSC. Here, long-wave irradiation, short-wave radiation (global radiation), wind direction and speed, amb temperature, precipitation and relative humidity were recorded at 1-minute inte starting in January 2021. Data gaps and time periods before January 2021 were filled data from two other climate measurement systems of the German weather service (D delivering open-source data within a range of 15 km from Bad Nauh (https://opendata.dwd.de/, accessed 27 January 2023).

Quasi-one-dimensional Simulations with FEFLOW ®
The simulations to investigate the influence of different groundwater level moisture content were performed using the FEFLOW® (v8.0) finite-element simula system from DHI Wasy GmbH [31]. A quasi-one-dimensional approach was adopte a detailed investigation at a specific site focusing on vertical water movement. For purpose, the models were created as one-meter-wide, two-dimensional vertical pl projections with a quadrangle element width of 1 m and height of 0.01 m. Water flow or out of the model is only possible through the surface or through the gravel aquifer schematic model setup is shown in detail in Figure 3. First, the model setup was validated by comparing the simulated moisture con with recorded values at the "Reference" and "mLSC" measuring points at 1.5 m, 2 and 2.9 m b.g.l. in time period 1 (b). Due to the expected significantly increased w extraction by the conifers at "Reference", the moisture sensors between 0.1 m and 0. b.g.l. were not considered. The value of the hydraulic head of the first time step wa as initial condition for all nodes.
At monitoring wells B1 to B5, the seasonal, spatial, and depth-dependent mois content differences over the course of different groundwater levels and surface influe were then determined by cyclic simulations in time period 2. For this purpose moisture contents of one simulated year at 1.5 m, 2.2 m, 2.9 m and 3.4 m b.g.l. reaching an annual equilibrium were considered. B1 to B5 are located outside the LSC are simply assumed to be the outer boundary of the LSC.  First, the model setup was validated by comparing the simulated moisture contents with recorded values at the "Reference" and "mLSC" measuring points at 1.5 m, 2.2 m and 2.9 m b.g.l. in time period 1 (b). Due to the expected significantly increased water extraction by the conifers at "Reference", the moisture sensors between 0.1 m and 0.75 m b.g.l. were not considered. The value of the hydraulic head of the first time step was set as initial condition for all nodes.
At monitoring wells B1 to B5, the seasonal, spatial, and depth-dependent moisture content differences over the course of different groundwater levels and surface influences were then determined by cyclic simulations in time period 2. For this purpose, the moisture contents of one simulated year at 1.5 m, 2.2 m, 2.9 m and 3.4 m b.g.l. after reaching an annual equilibrium were considered. B1 to B5 are located outside the LSC but are simply assumed to be the outer boundary of the LSC.

Material Properties and Model Layer Structure
The variably saturated states of the soils were simulated in FEFLOW ® using Richards' equation (mixed head-saturation) with the parametric van Genuchten-Mualem model [32][33][34]: L : parameter for pore tortuosity and pore connectivity.
The hydraulic material properties of the fine-grained layers were determined in advance using the HYPROP ® 2 measuring device from Meter Group, Inc. The measurement and evaluation are based on the simplified evaporation method [35][36][37]. Water retention and hydraulic conductivity curves, also based on the van Genuchten-Mualem model, were fitted by the automated adjustment of crucial parameters using the LABROS SoilView-Analysis software (v5.1.0.0), considering Equations (3) and (4). The van Genuchten-Mualem model was used to allow the direct transfer of the HYPROP ® results to FEFLOW ® .
A list of the soil samples used for the HYPROP ® tests and a graphical representation of the results are shown in Table A1 and Figure A1.
Basically, an alternation of silt loam and silty clay loam was assumed for simplification for the shallow fine-grained area based on soil profiles and grain size analyses from Rammler et al. [22].
The layered structure of the model for the "Reference" measuring point was determined using trial pit TP3. The corresponding borehole profiles were used for the models at monitoring wells B1 to B5. The soil composition of the measuring point "mLSC" is unknown and was derived from average values. Local variations, as well as the possible mixing of stratification due to excavation and backfilling in the LSC area, were neglected to ensure the comparability of the models. Furthermore, the layer boundaries might be less visible or disturbed in the drilled soil material. Bench gravel was divided into the less gravelly transition zone ("loam/gravel") and the pure gravel layer. Silty intermediate layers were neglected.
The layering for each model is summarized in Table 2. The material properties assigned to the layers based on the HYPROP ® results and literature data are listed in Table 3. The k 0 values determined with the Hyprop ® tests may differ from the actual saturated hydraulic conductivities [41][42][43]. In addition, the focus of the fitting process using LABROS SoilView-Analysis software was on the water retention curves. The porosities are based on the bulk densities of the soil samples, so deviations are possible here as well.
In the first step, the parameters k 0 , α and n of the loam/gravel layer were estimated with the FePEST (v8.0.2.3519) software by DHI Wasy GmbH using singular value decomposition based on the observed moisture contents at the "Reference" measuring point. The fitted material properties of this layer were then applied to the other models.

Aquifer and Surface Boundary Conditions
The recorded and manually measured groundwater levels at the respective monitoring wells were used as boundary conditions for the hydraulic head in the gravel aquifer. For the "Reference" measuring point, the data of monitoring well B6 were used, as it is located directly next to it. At the "mLSC" site, the exact groundwater level is not known, so average values were calculated from the data from B1 to B4. Due to the significantly different groundwater levels in the individual monitoring wells, their influence on the moisture content can be worked out by simulations. For the period from January to October 2020, the corresponding groundwater levels of 2021 were assumed.
On the surface of the models, the climatic influences were implemented as fluid flux boundary conditions ( Figure 3). These boundary conditions were derived from the basic water balance equation [39]: Due to the low slope, it was assumed that the subsurface discharge rate D is equal to groundwater recharge. For simplification, the potential evapotranspiration calculated with the FAO grass reference evapotranspiration model [44,45] based on the Penman-Monteith model [46] was assumed as the actual evapotranspiration.
The groundwater recharge rate at the test site is between 26 and 100 mm/a according to a rough graphical evaluation of groundwater recharge maps [47,48], with between 65% and 70% of precipitation evaporating in Hesse [49]. To obtain the corresponding values and variance of surface influences, two scenarios were considered, each assuming a certain percentage of precipitation as surface runoff, including additional evaporation of accumulated rainwater and plant interception (Table 4). By comparing the two scenarios, not only the influence of groundwater but also the surface influence on moisture content can be determined, and future developments due to climate change can be estimated.
No scenario distinction was made for the years 2020 and 2021, which were used in time period 1. For the dry year 2020, according to weather data, there was no groundwater recharge in the annual balance, and in the rainier and colder year 2021, after adjusting the surface runoff rate, there was an assumed groundwater recharge of 99 mm/a.

Derivation of Thermal Conductivities from Simulated Moisture Contents and Exemplary Simulations of Heat Extraction
In order to evaluate the influence of the different groundwater tables at the test site on the geothermal potential, thermal conductivities were derived from simulated moisture contents. For this purpose, two representative punctual thermal conductivity measurements from Rammler et al. [22] were used for the fine-grained soils for the investigated moisture ranges (θ > 0.15 cm 3 /cm 3 ) ( Table 5). Linear interpolation and extrapolation of the punctual measurements yielded the maximum possible and theoretical change in λ with changes in moisture content. However, the equations listed in Table 5 are not suitable for general thermal conductivity calculations, as they were only determined for a specific moisture content range.
The silty clay loam at the test site can have different thermal conductivities. For this study, a series of measurements with lower thermal conductivities were selected to obtain possible maximum differences from the silt loam layers.
For the loam/gravel layer, exemplary values from a master's thesis, in which the thermal conductivities of the gravelly soils of Bad Nauheim were investigated [50], were used (Table 5 and Figure A2b).
The measurement series are shown graphically in Figure A2. Table 5. Calculation equations of thermal conductivities (λ) using moisture content (θ) or saturation (S) for different depths and soil texture classes based on the thermal conductivity measurements shown in Figure A2. Also shown is the bulk density of the soil samples (ρ b ) during the measurement. For the thermal conductivity of the fine-grained soils, only the moisture range >0.15 cm 3 /cm 3 was considered. No grain size analysis is available for the silt loam soil sample, but it is assumed to be the same soil texture class. For the loam/gravel layer, λ is calculated based on saturation, since the bulk density might differ significantly from the value considered in the simulations. The potential impacts of the minimum and maximum groundwater-dependent or soil-texture-dependent thermal conductivities of the fine-grained soils calculated using the equations shown in Table 5 were illustrated by investigating the potential resulting heat extraction. Therefore, exemplary simulations were performed using the DELPHIN finitevolume simulation software (v6.1.5, https://www.bauklimatik-dresden.de/, accessed on 26 May 2023). With the DELPHIN software, the heat extraction of different geometric configurations of ground heat exchangers and LSCs can be examined in detail. Open-source climate data from the German weather service (https://opendata.dwd.de/, accessed on 27 January 2023) for moderate-climate regions in Germany were used for these simulations as boundary conditions. In addition, a homogeneous loamy subsurface with a constant thermal conductivity was assumed in each case. By adjusting the thermal conductivity of the soil according to the calculated minimum and maximum thermal conductivities, the resulting differences in heat extraction were determined.
The purpose of these simulations is to illustrate the noticeable variations in possible heat extraction associated with the variation in thermal conductivity. The results are therefore not presented as absolute values but as maximum percentage deviations in heat extraction that can be caused by different groundwater levels or the thermal properties of soils. The performed simulations do not correspond to an energy analysis of the specific LSC in Bad Nauheim.

Parameter Estimation and Model Validation at "Reference" Measuring Point
In the considered period from 7 November 2021 to 6 February 2023, no data on observed moisture contents at "Reference" were available on nine days (2%). The groundwater level fluctuated between 2.0 m and 3.6 m b.g.l. during this period, with a median of 3.1 m b.g.l.
The results of the simulations and the comparison with the observed moisture contents are shown in Figure 4. At 2.9 m b.g.l., the gravelly and loamy transition zone is present. Based on the material properties in Table 3, the simulation resulted in an RMSE value of 0.09 compared to the observed moisture contents (Figure 4c). By estimating the parameters k 0 , α and n with FePEST for this layer, the RMSE was significantly reduced to 0.01 due to a lower water retention capacity, regardless of the considered surface influences.
The adjusted material properties for the loam/gravel layer are shown in Table 6. Parameter estimation had no significant effect on the simulated moisture content at 2.2 m b.g.l. (Figure 4b). Low RMSE values of 0.02 and 0.03 resulted, with the simulated moisture contents mostly slightly higher than observed values. Due to the sharp rise in groundwater at the end of the simulation period (b), the silt loam was temporarily fully saturated at this depth. At 1.5 m b.g.l., the simulated moisture contents were significantly higher than the observed data, with RMSE values of 0.33 (Figure 4a). The median of the observed moisture contents was 0.07 cm 3 /cm 3 . The low values can be explained by the increased water extraction by the conifers in this area, which has an effect at least down to 1.5 m b.g.l. The resulting lower rate of groundwater recharge may also be the reason for the slightly higher simulated moisture contents at 2.2 m b.g.l. compared to observed values.
The low RMSE values of ≤0.03 at 2.2 m and 2.9 m b.g.l. in the range of measurement accuracy of the moisture sensors demonstrated the suitability of the quasi-one-dimensional model setup and the assigned material properties for the simulation of moisture contents, which, in this area, are mainly groundwater-influenced. Only minor or no differences were determined between scenarios A and B.
conductivities, the resulting differences in heat extraction were determined.
The purpose of these simulations is to illustrate the noticeable variations in possible heat extraction associated with the variation in thermal conductivity. The results are therefore not presented as absolute values but as maximum percentage deviations in heat extraction that can be caused by different groundwater levels or the thermal properties of soils. The performed simulations do not correspond to an energy analysis of the specific LSC in Bad Nauheim.

Parameter Estimation and Model Validation at "Reference" Measuring Point
In the considered period from November 07, 2021, to February 06, 2023, no data on observed moisture contents at "Reference" were available on nine days (2%). The groundwater level fluctuated between 2.0 m and 3.6 m b.g.l. during this period, with a median of 3.1 m b.g.l.
The results of the simulations and the comparison with the observed moisture contents are shown in Figure 4. At 2.9 m b.g.l., the gravelly and loamy transition zone is present. Based on the material properties in Table 3, the simulation resulted in an RMSE value of 0.09 compared to the observed moisture contents (Figure 4c). By estimating the parameters k0, α and n with FePEST for this layer, the RMSE was significantly reduced to 0.01 due to a lower water retention capacity, regardless of the considered surface influences.
The adjusted material properties for the loam/gravel layer are shown in Table 6.

Model Validation at Measuring Point "mLSC"
In the considered period from November 07, 2021, to February 06, 2023, no data on observed moisture contents were available at "mLSC" on 18 days (4%). The assumed groundwater table fluctuated between 3.7 m and 5.4 m b.g.l. during this period, and the median was 4.9 m b.g.l.
The results of the simulations and the comparison with the observed moisture contents are shown in Figure 5.
Between January and May 2022, the recorded moisture contents showed a rapid and significant decrease and, toward the end of this period, a rapid and significant increase at each depth. During this period, the soil temperatures also reached minimum values as low as −0.4 • C. It is therefore assumed that the soil was temporarily and partially frozen due to heat extraction by the LSC. The course of the moisture contents thus did not correspond to the undisturbed and natural course. The time delays in freezing and thawing in the comparison of moisture contents and soil temperature can be explained by the required withdrawal and supply of icing enthalpy. Judging from temperature trends, the soil partially refroze immediately after the end of the simulation period. Since the thermal effects of the LSC and thus possible soil freezing were not considered in the simulations, these specific periods were not used for validation.
RMSE values ≤ 0.04 were reached at all depths considered, with no differences between scenarios A and B. At 1.5 m b.g.l., a significant and rapid increase in moisture content was initially recorded just before freezing (Figure 5a). This was probably due to a measurement error caused by a previously insufficient connection of the sensor to the surrounding soil. Considering only the period after thawing, there were even smaller deviations, with an RMSE = 0.01. moisture contents, which, in this area, are mainly groundwater-influenced. Only minor or no differences were determined between scenarios A and B.

Model Validation at Measuring Point "mLSC"
In the considered period from November 07, 2021, to February 06, 2023, no data on observed moisture contents were available at "mLSC" on 18 days (4%). The assumed groundwater table fluctuated between 3.7 m and 5.4 m b.g.l. during this period, and the median was 4.9 m b.g.l.
The results of the simulations and the comparison with the observed moisture contents are shown in Figure 5. at "mLSC". Also shown is the root mean square error (RMSE) when comparing simulated and observed values. In addition, the soil temperature recorded by the temperature sensors is shown for each depth. The time period in which the soil was partially frozen and that was not considered in the RMSE calculation is indicated by the dashed lines. Partly due to only minor differences between scenarios A and B, the corresponding lines largely overlap.
Between January and May 2022, the recorded moisture contents showed a rapid and significant decrease and, toward the end of this period, a rapid and significant increase at each depth. During this period, the soil temperatures also reached minimum values as low as -0.4 °C. It is therefore assumed that the soil was temporarily and partially frozen at "mLSC". Also shown is the root mean square error (RMSE) when comparing simulated and observed values. In addition, the soil temperature recorded by the temperature sensors is shown for each depth. The time period in which the soil was partially frozen and that was not considered in the RMSE calculation is indicated by the dashed lines. Partly due to only minor differences between scenarios A and B, the corresponding lines largely overlap.
The largest deviations between simulated and observed moisture contents in the comparisons were observed at 2.9 m b.g.l. (Figure 5c), which is probably due to the unknown actual groundwater level. The significantly higher observed moisture content after thawing could be due to the appearance of cavities in the immediate vicinity of the sensor.
The low RMSE values ≤ 0.04, despite the groundwater levels being only assumed, were taken as evidence for the suitability of the quasi-one-dimensional model setup, including the assumed surface influences, and of the assigned material properties for the moisture content simulations.

Simulated Moisture Contents with Different Groundwater Levels at Monitoring Wells B1 to B5
The simulated site-and depth-dependent annual trends in moisture contents are shown in Figure 6a-d, exemplary of scenario A. Since the simulated values are based on an annually repeating cycle of groundwater level and surface influences as well as on an assumed soil structure, the actual moisture contents may deviate from them. However, the validation simulations at "Reference" and "mLSC" showed only minor deviations from the actually observed moisture contents. Basic trends, as well as the effects, of different groundwater conditions can therefore be derived from the simulated data. Since B1 to B5 are located at the outer edge of the LSC, and B1, in particular, is even farther away, the simulated moisture contents thus represent maximum and minimum values for the LSC area.
At B3 and B4, the groundwater level is temporarily higher than 2.9 m b.g.l. At B5, the highest groundwater table in time period 2 is 2.1 m b.g.l.
Based on the medians in both scenarios, moisture contents between 0.32 and 0.42 cm 3 /cm 3 were simulated for the fine-grained soils at the considered depths. The medians of the two scenarios differ by a maximum of 0.02 cm 3 /cm 3 . The sensors at "Reference" and "mLSC" recorded medians between 0.34 and 0.39 cm 3 /cm 3 , when periods of partially frozen soil are excluded. Thus, the observed values are within the range of moisture contents simulated at B1 to B5.
Moisture contents measured in a trial pit in the LSC area and moisture contents derived from gravimetric water contents measured in the boreholes B2 to B5 gave mean values between 0.31 and 0.34 cm 3 /cm 3 for fine-grained soils and are thus basically in a similar, slightly lower range of values [22]. Sampling for this evaluation was conducted in the summer and autumn of the dry year 2020. As can be seen in Figure 6, lower moisture contents generally occur during these seasons.  Based on the medians in both scenarios, moisture contents between 0.32 and 0.42 cm 3 /cm 3 were simulated for the fine-grained soils at the considered depths. The medians of the two scenarios differ by a maximum of 0.02 cm 3 /cm 3 . The sensors at "Reference" and "mLSC" recorded medians between 0.34 and 0.39 cm 3 /cm 3 , when periods of partially frozen soil are excluded. Thus, the observed values are within the range of moisture contents simulated at B1 to B5.
Moisture contents measured in a trial pit in the LSC area and moisture contents derived from gravimetric water contents measured in the boreholes B2 to B5 gave mean values between 0.31 and 0.34 cm 3 /cm 3 for fine-grained soils and are thus basically in a similar, slightly lower range of values [22]. Sampling for this evaluation was conducted in the summer and autumn of the dry year 2020. As can be seen in Figure 6, lower moisture contents generally occur during these seasons. At 3.4 m b.g.l., at B4 and B5, significantly lower moisture contents relative to the groundwater table occur temporarily due to the lower water retention capacity of the loam/gravel layer. In addition, different hydraulic properties of silt loam and silty clay loam are evident at a depth of 1.5 m b.g.l. While medians of 0.32 and 0.34 cm 3 /cm 3 were simulated in the silt loam at B1, the median values in the silty clay loam at other monitoring wells are significantly higher (up to 0.09 cm 3 /cm 3 ). Figure 6e-h show the differences between the simulated moisture contents of scenarios A and B. These generally decrease with depth, as the groundwater influence increases and the surface influence decreases. In addition, the differences are greatest at high groundwater distances, such as at B1 and B2. The maximum values occur at 1.5 m b.g.l., mainly in the rainy period at the beginning of the year, and then appear with a certain time delay as the depth increases. At 3.4 m b.g.l., at B4 and B5, significantly lower moisture contents relative to the groundwater table occur temporarily due to the lower water retention capacity of the loam/gravel layer. In addition, different hydraulic properties of silt loam and silty clay loam are evident at a depth of 1.5 m b.g.l. While medians of 0.32 and 0.34 cm 3 /cm 3 were simulated in the silt loam at B1, the median values in the silty clay loam at other monitoring wells are significantly higher (up to 0.09 cm 3 /cm 3 ). Figure 6e-h show the differences between the simulated moisture contents of scenarios A and B. These generally decrease with depth, as the groundwater influence increases and the surface influence decreases. In addition, the differences are greatest at high groundwater distances, such as at B1 and B2. The maximum values occur at 1.5 m b.g.l., mainly in the rainy period at the beginning of the year, and then appear with a certain time delay as the depth increases.

Influence of Spatial Differences in Groundwater Distance on Moisture Contents
The comparison of the simulated moisture contents with the distance to the groundwater level at the respective depths shows a clear correlation in the form of a power function with coefficients of determination (R 2 ) between 0.86 and 0.99 (Figure 7).
Moisture content decreases with increasing distance to the groundwater level, depending on the soil texture class and surface influences. Considering medians (θ med ), differences of up to 0.08 cm 3 /cm 3 occur within one soil texture class (here: silt loam).
No differences between scenarios A and B are observed at very low groundwater distances. Here, the moisture content is mainly determined by capillary rise and groundwaterrelated temporary full saturation. With increasing distance to groundwater, increasing differences in θ med between the scenarios of up to 0.02 cm 3 /cm 3 occur with generally higher moisture contents in scenario A, which is due to the higher groundwater recharge rate considered in scenario A. In addition, the decrease in moisture content with increasing groundwater distance in scenario B is therefore steeper. The comparison of the simulated moisture contents with the distance to the groundwater level at the respective depths shows a clear correlation in the form of a power function with coefficients of determination (R 2 ) between 0.86 and 0.99 (Figure 7). Moisture content decreases with increasing distance to the groundwater level, depending on the soil texture class and surface influences. Considering medians (θmed), differences of up to 0.08 cm 3 /cm 3 occur within one soil texture class (here: silt loam).
No differences between scenarios A and B are observed at very low groundwater distances. Here, the moisture content is mainly determined by capillary rise and groundwater-related temporary full saturation. With increasing distance to groundwater, increasing differences in θmed between the scenarios of up to 0.02 cm 3 /cm 3 occur with generally higher moisture contents in scenario A, which is due to the higher groundwater recharge rate considered in scenario A. In addition, the decrease in moisture content with increasing groundwater distance in scenario B is therefore steeper.
For the same depth and the same soil texture class, differences in moisture content of up to 0.08 cm 3 /cm 3 result from the spatial variation in groundwater levels. The largest differences occur at 2.9 m b.g.l. in the silt loam layer.

Influence of Groundwater on the Seasonal Variation in Moisture Content
The comparison of the maximum annual range of moisture contents (Δθ) at a given depth and location with θmed for different partially saturated soil texture classes and depth ranges shows a clear linear correlation, with Pearson's correlation coefficients (r) between 0.87 and 0.99 (Figure 8). For the same depth and the same soil texture class, differences in moisture content of up to 0.08 cm 3 /cm 3 result from the spatial variation in groundwater levels. The largest differences occur at 2.9 m b.g.l. in the silt loam layer.

Influence of Groundwater on the Seasonal Variation in Moisture Content
The comparison of the maximum annual range of moisture contents (∆θ) at a given depth and location with θ med for different partially saturated soil texture classes and depth ranges shows a clear linear correlation, with Pearson's correlation coefficients (r) between 0.87 and 0.99 (Figure 8). Depending on the depth and the soil texture class, Δθ increases with increasing median moisture content and thus with decreasing distance to groundwater, additionally depending on the groundwater recharge rate. Maximum values of 0.08 cm 3 /cm 3 are reached. In general, higher fluctuations are observed in scenario A due to the higher recharge rate from precipitation. At groundwater distances < 2 m from the considered depth, no differences between the two scenarios can be observed due to the groundwater influence.
Due to the surface influence of precipitation, high values of Δθ generally occur at 1.5 m b.g.l. At 2.2 m b.g.l., significantly higher fluctuations occur at B5, since the groundwater rises to 2.1 m b. g. l. at times, and thus fully saturated conditions also prevail. Depending on the depth and the soil texture class, ∆θ increases with increasing median moisture content and thus with decreasing distance to groundwater, additionally depending on the groundwater recharge rate. Maximum values of 0.08 cm 3 /cm 3 are reached. In general, higher fluctuations are observed in scenario A due to the higher recharge rate from precipitation. At groundwater distances < 2 m from the considered depth, no differences between the two scenarios can be observed due to the groundwater influence.
Due to the surface influence of precipitation, high values of ∆θ generally occur at 1.5 m b.g.l. At 2.2 m b.g.l., significantly higher fluctuations occur at B5, since the groundwater rises to 2.1 m b.g.l. at times, and thus fully saturated conditions also prevail. At 2.9 m and 3.4 m b.g.l., fully saturated conditions exist intermittently at B3, B4 and B5 and additionally at 2.2 m b.g.l. at B5. Therefore, high ∆θ and θ med values for silt loam are generally observed here. However, the range of variation decreases with increasing moisture content and thus with decreasing groundwater distance (r = −0.84). Since high moisture contents are already present due to the very low groundwater distance, the difference between the moisture content in the temporarily partially saturated conditions and the fully saturated moisture content is lower.

Derived Thermal Conductivities from Simulated Moisture Contents at Monitoring Wells B1 to B5
The derived site-and depth-dependent annual trends in thermal conductivity are shown in Figure 9a-d. Since the derived values are based on an annually repeating cycle of simulated moisture contents and on exemplary laboratory measurements of λ, the actual thermal conductivities may deviate from them. In addition, temporary icing caused by the LSC, as can be seen in the observed moisture contents in Figure 5, was not considered. Frozen soils generally have a significantly increased thermal conductivity [51]. (e) (f) (g) (h) Figure 9. (a-d) Thermal conductivities derived from the simulated moisture contents with the equations from Table 5 for scenario A and (e-h) the difference in the derived thermal conductivities for scenario A (λA) and scenario B (λB) at monitoring wells B1 to B5 and their specific groundwater levels and soil structures at the depths 1.  Table 5 for scenario A and (e-h) the difference in the derived thermal conductivities for scenario A (λ A ) and scenario B (λ B ) at monitoring wells B1 to B5 and their specific groundwater levels and soil structures at the depths 1.5 m, 2.2 m, 2.9 m and 3.4 m b.g.l. Only fine-grained silty layers are considered for (a-d). Time = 0 corresponds to January 01.
Medians of λ between 1.52 W/(m·K) at B1 and 1.68 W/(m·K) at B5 were derived for the silt loam layers at the considered depths. Due to the deviating thermal properties of the silty clay loam, the medians in this layer at 1.5 m b.g.l. range between 1.32 W/(m·K) and 1.35 W/(m·K) (Figure 9a). The derived values largely agree with the mean thermal conductivities determined by Rammler et al. [22].
Looking at the medians, differences of up to 0.21 W/(m·K) at 1.5 m b.g.l. result, which can mainly be attributed to the different thermal properties of silt loam and silty clay loam. Using the maximum differences at this depth with 1.34 W/(m·K) at B2 and 1.55 W/(m·K) at B1 in scenario A for the exemplary DELPHIN simulations, the different thermal properties of the fine-grained soils result in maximum spatial variations in theoretical heat extraction of 13% to 14% at this depth. Figure 9e-h show the differences between the derived values for scenarios A and B. As with the moisture content, the differences generally decrease with depth, as the groundwater influence increases and the precipitation influence decreases. Due to the comparatively high groundwater distances at B1 and B2, the differences are greatest there, with values of up to 0.05 W/(m·K). The maximum values occur mainly in the rainy period and then appear with a certain time delay as the depth increases. In the area with high groundwater levels (B3 to B5), there are no or only minor differences, particularly at 2.9 m and 3.4 m b.g.l.
The high ∆θ values of the loam/gravel layer also cause the derived thermal conductivity to vary by up to 0.92 W/(m·K) at 3.4 m b.g.l. (Figure 10). With a maximum thermal conductivity of 2.36 W/(m·K), higher values are obtained in this layer than in the fine-grained soils at the same depth by up to 0.82 W/(m·K) when comparing medians. For this layer, there is no significant difference between scenarios A and B due to the low groundwater distance and long-time full saturation.  Table 5 for (a) scenario A and (b) scenario B at monitoring wells B1 to B5 and their specific groundwater levels and soil structures at 3.4 m b.g.l. considering fine-grained silts and the loam/gravel layer. Time = 0 corresponds to January 01.

Influence of Spatial Differences in Groundwater Distance on the very Shallow Geothermal Potential (vSGP)
The comparison of the derived thermal conductivities with the distances to groundwater at the respective depths shows a clear correlation in the form of a power function with coefficients of determination (R 2 ) between 0.85 and 0.99 ( Figure 11).   Table 5 for (a) scenario A and (b) scenario B at monitoring wells B1 to B5 and their specific groundwater levels and soil structures at 3.4 m b.g.l. considering fine-grained silts and the loam/gravel layer. Time = 0 corresponds to January 01.

Influence of Spatial Differences in Groundwater Distance on the Very Shallow Geothermal Potential (vSGP)
The comparison of the derived thermal conductivities with the distances to groundwater at the respective depths shows a clear correlation in the form of a power function with coefficients of determination (R 2 ) between 0.85 and 0.99 ( Figure 11).
As the groundwater distance increases, λ med decreases due to the moisture content also decreasing (Figure 7). Looking at medians, the thermal conductivity of the silt loam soil texture class shows differences of up to 0.16 W/(m·K). As the distance to groundwater increases, the difference between scenarios A and B increases to 0.03 W/(m·K) when comparing median values, as the influence of groundwater on the moisture content decreases and that of precipitation increases. The resulting higher moisture contents in scenario A lead to higher thermal conductivities. The lower groundwater recharge rate also leads to a steeper decrease in λ med with increasing groundwater distance in scenario B. However, at low groundwater distances, the very shallow geothermal potential in the form of λ is mainly influenced by capillary rise and temporary full saturation.
For the same depth and the same soil texture class, differences in thermal conductivity of up to 0.14 W/(m·K) result from the spatial variation in groundwater levels, considering medians. The largest differences occur at 2.9 m b.g.l. in the silt loam layer. Using this maximum spread of thermal conductivity for the exemplary DELPHIN simulations, with 1.54 W/(m·K) at B1 and 1.68 W/(m·K) at B5 in scenario B, this results in a maximum spatial variation in heat extraction of 7% to 9%.
(a) (b) Figure 10. Thermal conductivities derived from the simulated moisture contents with the equations from Table 5 for (a) scenario A and (b) scenario B at monitoring wells B1 to B5 and their specific groundwater levels and soil structures at 3.4 m b.g.l. considering fine-grained silts and the loam/gravel layer. Time = 0 corresponds to January 01.

Influence of Spatial Differences in Groundwater Distance on the very Shallow Geothermal Potential (vSGP)
The comparison of the derived thermal conductivities with the distances to groundwater at the respective depths shows a clear correlation in the form of a power function with coefficients of determination (R 2 ) between 0.85 and 0.99 ( Figure 11). As the groundwater distance increases, λmed decreases due to the moisture content also decreasing (Figure 7). Looking at medians, the thermal conductivity of the silt loam soil texture class shows differences of up to 0.16 W/(m•K). As the distance to groundwater increases, the difference between scenarios A and B increases to 0.03 W/(m•K) when comparing median values, as the influence of groundwater on the moisture content decreases and that of precipitation increases. The resulting higher moisture contents in scenario A lead to higher thermal conductivities. The lower groundwater recharge rate also leads to a steeper decrease in λmed with increasing groundwater distance in scenario B. However, at low groundwater distances, the very shallow geothermal potential in the form of λ is mainly influenced by capillary rise and temporary full saturation.
For the same depth and the same soil texture class, differences in thermal conductivity of up to 0.14 W/(m•K) result from the spatial variation in groundwater levels, considering medians. The largest differences occur at 2.9 m b.g.l. in the silt loam layer. Using this maximum spread of thermal conductivity for the exemplary DELPHIN

Influence of Groundwater on the Seasonal Variation in the Very Shallow Geothermal Potential (vSGP)
A comparison of the maximum annual range of thermal conductivity (∆λ) at a given depth and location with λ med for different partially saturated soil texture classes and depths shows a clear linear correlation, with Pearson's correlation coefficients (r) between 0.87 and 1.00 ( Figure 12).

Influence of Groundwater on the Seasonal Variation in the very Shallow Geothermal Potential (vSGP)
A comparison of the maximum annual range of thermal conductivity (Δλ) at a given depth and location with λmed for different partially saturated soil texture classes and depths shows a clear linear correlation, with Pearson's correlation coefficients (r) between 0.87 and 1.00 ( Figure 12). Depending on the depth and soil texture class, Δλ increases with increasing moisture content and thermal conductivity and thus with decreasing distance to groundwater with a maximum of 0.14 W/(m•K). As with the moisture content, the variations are smaller for scenario B due to the lower recharge rate from precipitation. There is no difference at low groundwater distances, where the influence of the groundwater dominates.
Due to the groundwater-induced high moisture content and temporary full saturation, there is also a decrease in Δλ with increasing λmed for the affected soils (r = -0.85).
The maximum Δλ value of 0.14 W/(m•K) occurs, for example, at B5 at a depth of 2.2 m b.g.l. Here, the thermal conductivity varies between 1.57 W/(m•K) and 1.71 W/(m•K), mainly related to the annual groundwater fluctuation. Based on these values, the exemplary DELPHIN simulations result in a maximum annual variation in theoretical heat extraction of 7% to 8%. Depending on the depth and soil texture class, ∆λ increases with increasing moisture content and thermal conductivity and thus with decreasing distance to groundwater with a maximum of 0.14 W/(m·K). As with the moisture content, the variations are smaller for scenario B due to the lower recharge rate from precipitation. There is no difference at low groundwater distances, where the influence of the groundwater dominates.

Conclusions
Due to the groundwater-induced high moisture content and temporary full saturation, there is also a decrease in ∆λ with increasing λ med for the affected soils (r = −0.85).
The maximum ∆λ value of 0.14 W/(m·K) occurs, for example, at B5 at a depth of 2.2 m b.g.l. Here, the thermal conductivity varies between 1.57 W/(m·K) and 1.71 W/(m·K), mainly related to the annual groundwater fluctuation. Based on these values, the exemplary DELPHIN simulations result in a maximum annual variation in theoretical heat extraction of 7% to 8%.

Conclusions
The validation of the simulations based on the observed moisture contents at the measuring points "Reference" and "mLSC" showed the general suitability of the quasione-dimensional model setup and of the assigned material properties to simulate actual moisture contents. The methodology can thus also be transferred to other applications or issues concerning vertical water movement in a soil column. The extent to which scenario A or B corresponds more to the actual surface influences could not be conclusively assessed due to the small differences between the simulated moisture contents.
The simulations at monitoring wells B1 to B5 resulted in the following influencing factors and effects on moisture content and on the very shallow geothermal potential in the form of the thermal conductivity in the area around the LSC (Table 7).

∆ max (λ) [W/(m·K)]
Groundwater Spatial difference at the same depth for the same fine-grained soil texture class 0.07 0.14 Annual variation in fine-grained soils 0.08 0.14 Annual variation in loam/gravel layer at 3. According to exemplary simulations, groundwater-related spatial and seasonal variations in thermal conductivity can cause maximum changes in heat extraction of 7% to 9%. The determined differences represent maximum values due to the location of the monitoring wells at the outer boundary of the LSC.
Due to the different soil texture classes of fine-grained soils and their thermal properties, variations in heat extraction of 13% to 14% are possible. Even larger grain-sizedistribution-related differences in thermal conductivity are particularly evident with respect to the loam/gravel layer, which is present in the northeast of the LSC from 3 m b.g.l. downward. Particularly high thermal conductivities are also expected in the purer gravel [50,52], which is only immediately deeper. Due to the effective thermal conductivity, which is influenced by the groundwater flow velocity [53], the bench gravel as a whole can thus have a significant impact on the thermal regeneration capability of the subsurface, depending on the location at the LSC.
Accordingly, a sufficiently detailed subsurface exploration is required in advance for the design of an efficient LSC. In addition to determining the spatial and seasonal groundwater conditions, the distribution of different soil textures and specific thermal properties should be a primary focus. In practice, this can be achieved primarily through a higher number of subsurface outcrops corresponding to the large area of an LSC. However, this must also be determined site-specifically based on various factors, such as the slope of the terrain or on the basis of hydrogeological, geological and pedological maps. A sufficiently detailed subsurface exploration can potentially lead to smaller and more efficient collector and distribution systems, thus increasing overall profitability.
The investigation has also shown that as groundwater recharge from precipitation decreases because of climate change [53], differences in the very shallow geothermal potential between areas at the LSC with high (east) and low (west) groundwater tables will Table A1. List of soil samples for which hydraulic properties were determined using the HYPROP ® measuring device for each layer or soil texture class. The sampling sites are documented in Rammler et al. [22]. In addition, further tests were conducted that are not listed here because they were not used for this study. ρ b_ini corresponds to the bulk density before saturation. Appendix A Table A1. List of soil samples for which hydraulic properties were determined using the HYPROP® measuring device for each layer or soil texture class. The sampling sites are documented in Rammler et al. [22]. In addition, further tests were conducted that are not listed here because they were not used for this study. ρb_ini corresponds to the bulk density before saturation.  Figure A1. Water retention curves based on the Genuchten-Mualem model with soil moisture tension (pF) and moisture content (θ) fitted to the evaluated values of the soil samples from Table  A1 for (a) topsoil, (b) silt loam and (c) silty clay loam with LABROS SoilView-Analysis. The resulting parameters are listed in Table 3.  Table A1 for (a) topsoil, (b) silt loam and (c) silty clay loam with LABROS SoilView-Analysis. The resulting parameters are listed in Table 3.

Layer/Soil
(a) (b) (c) Figure A1. Water retention curves based on the Genuchten-Mualem model with soil moisture tension (pF) and moisture content (θ) fitted to the evaluated values of the soil samples from Table  A1 for (a) topsoil, (b) silt loam and (c) silty clay loam with LABROS SoilView-Analysis. The resulting parameters are listed in Table 3.  [22] on which the derivation of the thermal conductivity (λ) was based for (a) fine-grained soils as a function of moisture content (θ) and (b) the loam/gravel layer as a function of saturation (S) ( Table 5). The series of measurements selected were those of sample TP6 (1.6 m) for silty clay loam and those of sample TP2 (2.4 m) for silt loam. The sampling sites are documented in Rammler et al. [22]. Also shown in (a) are the trend lines of linear interpolation and extrapolation for θ > 0.15 cm 3 /cm 3 , as described in Table 5. For the loam/gravel layer, a series of measurements on the corresponding layer at B1 from a master's thesis [50] was used. The trend line based on a second-degree polynomial function is also shown. The coefficient of determination (R 2 ) is given for all trend lines.