Quantifying the contribution of heat recharge from confining layers to geothermal resources

Geothermal operations are expanding and increasingly contributing to the current energy supply. Assessing the long-term operable lifetime of these projects is complicated as the reservoirs they produce from are often deep and subsurface properties are uncertain and spatially variable. The minimum lifetime of a geothermal project usually considers the heat in place in a geothermal reservoir, but not the heat flow from the confining layers into the reservoir during operation. For the economic feasibility and optimal design of a geothermal project it is the key to capture this process, as this allows more accurate prediction of the long-term extraction temperature. Previous studies are not conclusive on the contribution of vertical recharge from the confining layers. This research evaluates the contribution of recharge to the heat production for geothermal projects. In a simulation study with an idealised, homogeneous geothermal system, the reservoir thickness, well placement and production rate are varied to investigate their respective influence on thermal recharge. The results show that the recharge from the vertical confining layers may contribute considerably to the total energy output of a geothermal well. Due to recharge, geothermal systems produce more heat than the heat in place and under specific conditions the produced heat can be more than five times as large. The largest contribution of the recharge of vertically confining layers can be expected under conditions of thin reservoirs and a long interaction time between the injection water and confining layers (i.e. large well spacing, low production rate). The insights of this study may help to optimise the design and operation of individual geothermal projects for optimal utilisation of a geothermal reservoir for energy supply.


Introduction
The Earths' interior contains enough heat to supply a considerable part of the worlds' energy need (Dickson and Fanelli, 2003;Rybach, 2015).Along the edges of continental plates and in volcanic areas, geothermal heat is available through high temperature and at shallow depths (Stefánsson, 1998;Carranza et al., 2008).However, in most regions with urban development where demand for geothermal heat exists, the thick crust reduces the heat flow from the Earths' mantle to the surface.As a result the possibility to utilise geothermal heat is challenging in such regions (Mijnlieff, 2020;Willems and Nick, 2019).
In spite of the challenges associated with the limited amount of heat being available, geothermal operations are expanding and contributing to the energy supply of these regions.For these projects, it is difficult to estimate the efficiency of their heat production, as the reservoirs they produce from are deep and subsurface properties are uncertain and spatially variable.With the uncertainties in the long-term extraction temperature of the geothermal well, the operable lifetime of geothermal projects becomes uncertain.Moreover, the demand for heat is mostly concentrated in urban or greenhouse areas where multiple users have time-varying needs for energy, which further complicates the estimation of the lifetime of a single geothermal project (Willems et al., 2017a,b).In practice, the minimum lifetime of a geothermal project is usually determined by taking into account the heat in place in a reservoir at the start of operation, together with the well placement and production rate (Willems and Nick, 2019;Daniilidis et al., 2017;Lopez et al., 2010).These calculations do not always take into account the heat flow from the confining layers that recharges the geothermal reservoir during exploitation.This heat flow affects the production of heat during and beyond, the lifetime of the geothermal project.As a worst-case approximation for the economic life-time of a single geothermal project this may be an acceptable way to go.But also after the economic life-time of a geothermal project still heat is needed, as cities will not move during that time.Therefore, for the economic feasibility, optimal design and long term operation of a geothermal project it is the key to capture heat recharge, as this allows for a more accurate prediction of the development of the long-term extraction temperature.Previous research explored which factors affect the long-term extraction temperature of geothermal projects, but the results on the contribution of vertical recharge from the confining layers are not conclusive or consistent.Also some studies ignore the role of recharge from vertically confining layers (Axelsson et al., 2005;Axelsson, 2012;Kong et al., 2017;Bauer et al., 2019;Crooijmans et al., 2016;Lopez et al., 2010).While others emphasise the role of recharge from the confining layers (Poulsen et al., 2015;Daniilidis and Herber, 2017;Randolph et al., 2011;Wang et al., 2020;Shetty et al., 2018).Not all of these studies distinguish between recharge from the vertical confining layers and lateral heat flow from outside the reservoir area.Given the potential interference between neighbouring geothermal projects in the same reservoir, it is important to differentiate between laterally and vertically recharged heat.Building on the previous work on thermal recharge, the goal of this research is to (A) evaluate the contribution of horizontal and vertical recharge to the heat produced by a geothermal system and (B) identify under which conditions this contribution is largest.
For this study, we use an idealised, homogeneous geothermal system in which reservoir thickness, well placement and production rate are varied to investigate their respective influence on thermal recharge.Ultimately, the results of this study may help developers and operators to optimise the contribution of recharge in their geothermal projects.

Fluid flow
The reservoir water is used as a carrier for the heat by extracting it with a tube well, taking out the heat via a heat exchanger and simultaneous re-injection of the reservoir brine in a second well.The groundwater flow through the reservoir matrix is described by the groundwater flow equation (single porosity, single phase and incompressible) (Langevin et al., 2007): = ρS s,0 ∂h 0 ∂t s is a source or sink of fluid with density ρ s .

Heat transfer processes
Due to the temperature difference between the injected water and the reservoir ambient temperature, heat transfer in the reservoir will take place.Heat transfer occurs in the form of convection (movement of the water) and conduction (initial heat of the reservoir rock).
The following equation describes heat transport in a porous medium (Langevin et al., 2007): The nomenclature provides each parameter of Eq. (2.2).The main terms in the equation are: : thermal retardation.This term reflects the delay of the cold front relative to the injected water, due to heat storage in the reservoirs' matrix (Bloemendal et al., 2018).
: conduction of heat.The conductive heat flux is directly proportional to a temperature gradient, expressed by Fourier's Law (Ferrell and Stahel, 2000): (2.3) this law states that conductive heat flows from a high temperature to a low temperature, and that the amount of heat flow is dependent on the conductivity of the medium through which it is travelling.
Thermal diffusivity is used to calculate the conductive heat transport with a diffusion equation.Thermal diffusivity describes the rate of temperature spread through a material (Ferrell and Stahel, 2000).Thermal diffusivity is determined by a combination of thermal conductivity and volumetric heat capacity.When more heat is flowing into a unit volume than flowing out of it, the temperature will rise.The rate of the heat flow is determined by the thermal conductivity.The temperature increase of the material is dependent on the heat capacity (Chekhonin et al., 2012).The thermal diffusivity term can be calculated using the following equation (Langevin et al., 2007): fluid heat capacity [J/kg/ (2.4) • α q θ : dispersion.Local variations in flow velocity result in a smooth transition between injected and ambient temperature.
• qT: advection, the heat transport in porous media by the flow of ground water.• q ′ s T s : heat transport to or from a sink/source in the model domain.

Thermal recharge
A geothermal reservoir is in open contact with the rest of the reservoir and its confining layers.When the reservoir cools down due to the injection of cold water, a temperature difference between the confining layers and the reservoir emerges.This results in a heat flow from the confining layers towards the reservoir.Due to the low permeability of the layers vertically confining the reservoir, pore fluid flow is minimal and the transfer of heat is based on conduction.This phenomenon, which is central in this study, is referred to as vertical thermal recharge.
In addition to the vertical recharge, heat will also flow laterally from within the reservoir to the cooled-down area around the injector well.This heat flow is driven by both conduction and free convection and called lateral thermal recharge.Free convection occurs when permeability is high.

Model setup
In this study, we use SEAWAT to simulate the production from a geothermal reservoir.It allows for incorporating all the required processes in an open-source flexible environment and it has the ability to incorporate various subsurface conditions in the simulation and monitoring of the heat flow within and into a geothermal reservoir (Bakker et al., 2016).SEAWAT is a well-established model and it is used in various similar studies, e.g.(Zeghici et al., 2015;Hecht-Méndez et al., 2010;Vandenbohede et al., 2014).SEAWAT is a finite-difference computer code designed to simulate coupled three-dimensional, variable-density groundwater flow and heat transport (Guo et al., 2002;Thorne et al., 2006).SEAWAT uses Eq. (2.2) to determine the heat transport in the model.The default SEAWAT configuration assumes a linear temperature-density dependency (Langevin et al., 2008), which can be done for systems with small temperature changes.However, in this study the temperature changes are too large, therefore, we include a polynomial approximation of the temperature-density dependency, as was done by van Lopik et al. (2016), Rocchi (2019) and Marif (2019).
Model domain and characteristics.In our study, we represent a homogeneous geothermal system with a box model that simulates a reservoir at a depth of 2500 m.The reservoir extends laterally over the entire model, and has a size of 6.2 km × 6.2 km.With exploratory simulations it was identified that extending the model domain beyond 6.2 × 6.2 km did not affect the simulation results, indicating that at this domain size the model boundaries do not affect the result anymore.An 810 m thick over and under burden confine the reservoir.
The properties of the reservoir material are set to match the properties of a sandstone as we observe it in the geology of The Netherlands, the properties of the confining layers material are set to match the properties of a shale.Following representative Dutch field data, the applied parameter values (International Association for the Properties of Water and Steam, 1994) are provided in Table 1.Since salinity is constant throughout the model domain, the volumetric heat capacity only varies due to temperature change, as a result of density changes (<2%).Since these are very small changes, the volumetric heat capacities are constant in the model.
The outer cells of the model are all set with a constant temperature and head boundary.A constant head/temperature boundary provides an inexhaustible source or sink of water/heat to and from the model domain.The used heat flow from the centre of the Earth towards the subsurface, is 65 mW/m 2 (Davies and Davies, 2010) by assigning an initial temperature in the model domain at the boundary, the applied gradient is provided in Table 1.In all experiments, we assume a single geothermal doublet that consists of vertical wells with a perforated screen along the entire thickness of the reservoir.The injection temperature is constant throughout the entire production period at a temperature of 35 • C.
Spatial and temporal discretisation.To provide stable and accurate results with the SEAWAT groundwater flow and heat-transport model, spatial and temporal discretisation is chosen such to minimise numerical errors, while at the same time maintain low computational effort.
The model consists of 100 by 100 cells in the horizontal x-and y-direction with a size of 40 m by 40 m in the direct vicinity of the doublet (600 m around the wells).To minimise the effect of the boundary conditions on the flow in the reservoir, and at the same time minimise the computational effort, the cells beyond the area of interest have a size of 100 m by 100 m.In the vertical direction the model consists of a reservoir confined by an over-and an under burden.The model layers in the reservoir have a thickness of 10 m, the number of model layers changes accordingly with the reservoir thickness scenario.The vertical confining layers are always the same: both are made up of 18 model layers, of which the most outer one and the 10 model layers adjacent to the reservoir have a thickness of 10 m, the remaining 7 model layers have a thickness of 100 m.
Spatial discretisation must be consistent with the dispersivity parameters (Reyes et al., 2013).The injection and production well are placed at the same y-coordinate of the model, therefore the main flow direction is in the y-direction.Longitudinal dispersivity, α L , is used to represent the local variations in the velocity field of a groundwater solute in the direction of fluid flow (Shlomo, 2006).Dispersion perpendicular to the flow direction is represented by horizontal and vertical transverse dispersivity α Th and α Tv (Delgado, 2007).The longitudinal dispersivity in this model is 10 m.The transverse dispersivity tensors are represented by the ratio of the transverse dispersivity to the longitudinal dispersivity.For the horizontal and vertical transverse dispersivity this ratio is respectively 0.1 and 0.01 (Vandenbohede et al., 2014;Zeghici et al., 2015).This means that the longitudinal dispersion and mixing dominates over the transvere dispersion and mixing, ensuring that vertical heat flow is not influenced by any dispersion related effects.The horizontal cell size and the longitudinal dispersion value result in a Peclet number Pe of 4 for this model (Poulsen et al., 2015).Literature argues that only a Peclet number smaller than 2 guarantees numerical stability (Reyes et al., 2013), and with a Peclet number smaller than 4 it can be expected (Poulsen et al., 2015).In the heat transport equation solver, SEAWAT reduces time steps by itself in order to meet the condition set.
Validation.To check whether the model with the given grid size and Pe = 4 is numerical stable, we conduct a set of experiments in which the cell size in the reservoir and confining layers is decreased step by step until no more changes occurred in the subsequent simulation results.These results indicated that the model is numerically stable with the used cell sizes, and has limited numerical dispersion.SEAWAT automatically adjusts its time step to meet the Courant condition which was set to 0.8.Furthermore, our model set-up is the same as the work carried out by (Poulsen et al., 2015).The results of both models showed a good match, which suggests that our model captures the relevant processes and produces realistic results.This allows us to compare different scenarios with the goal to identify the absolute and relative contribution of recharge on the production temperature and heat content of the reservoir.Due to the lack of historical field data, we could not calibrate the model to real data.

Experiments
The starting point of this study is the simulation of two cases, with and without the confining layer activated in the model.The goal of these simulations is to show to what extent the vertical recharge affects performance of geothermal projects.These cases have the following properties: For both cases well spacing (L) is 1200 m, and injection temperature (T inj ) is 35 • C. To compare the effect of a change in reservoir thickness, the production rate increases proportionally with reservoir thickness (H).
In a second set of experiments the reservoir thickness, production rate and well spacing are varied one-by-one to assess their influence on the recharge: • Reservoir thickness: H = 50, 100, 200 m • Production rate: Q = 150, 300, 600 m 3 /h (Q = 0.042, 0.083, 0.17 m 3 /s) • Well spacing: L = 1000, 1200, 1500 m The two properties that are not assessed remain stable at H = 200 m, Q = 300 m 3 /h and L = 1200 m.The characteristics of the scenarios are presented together with the results in one overview table in the results section (Table 2).These experiments all included vertical confining layers and allow for vertical recharge.The model simulates the heat production from a geothermal well for 300 years, with a constant flow rate and re-injection temperature.The 300-year period was chosen to allow for assessment of the long-term effect of recharge, beyond the economic life time of geothermal projects.It also allowed for straightforward comparison with the results of Poulsen et al. (2015) for model validation.

Assessment framework
During the simulation period the temperature in the production well and the total energy content of the reservoir and confining layers are monitored.Differences in the production temperature between different experiments illustrate the performance of the geothermal project.The time-integrated extracted energy from the production well [J], E p (2.5) Bruijn et al. is calculated as follows: (2.5) with t the time in years with starting point t 0 , ΔT av the average difference between the production temperature and the injection temperature, c f the fluid volumetric heat capacity [J/m 3 /K] and Q the production rate [m 3 /h].In all the results, the production temperature is in fact the down-hole temperature in the model.To allow for discrimination between vertical and lateral recharge and to assess the heat production relative to the heat in place, a representative reservoir volume (V [m 3 ]) is defined following the commonly applied method in the Netherlands for demarcation of exploration licences (Mijnlieff and Van Wees, 2009;de Zwart et al., 2012;Mijnlieff, 2014).The representative reservoir volume is defined based on the well spacing: in which H is the reservoir thickness [m], and L is the well spacing [m] (Fig. 1).
The initial energy content of the representative volume E rep [J], is determined by: The amount of thermal energy that is recharged during the simulation period towards the representative volume E r , is calculated using the method introduced by Poulsen et al. (2015), which states that When ΔE rep,t is the time-integrated change in energy content of the representative volume at time t, E p the time-integrated extracted energy from the production well, then E r is the amount of thermal energy that is recharged during production towards the representative volume:   representing the confining layers.As the total recharge consists of vertical and lateral recharge, the laterally recharged heat, E rl follows from the difference between the total E r and vertical recharge E r,v : (2.9) The initial representative reservoir heat content, the produced and recharged heat are used to determine dimensionless production ratios to assess the relative effect of the recharge.
Firstly, the produced energy ratio ϵ p is defined as: Secondly, the heat vertically recharged ratio ϵ rv is defined as: Lastly, the recharged heat into the representative volume ratio, both vertically and laterally, ϵ r is defined as: (2.12)

Influence of the confining layers
We investigate the effect of the confining layers for the models with reservoir thickness of 50 and 200 m, and a proportional production rate of respectively Q = 150 and 600 m 3 /h.Fig. 2 shows the production temperature as a function of time.The result of the simulation with reservoir thickness of 50 m shows significantly less temperature decline, while the 200-m case gives results similar to those in the no-confininglayers case.If heat flow in the confining layers is disabled, as indicated with the dashed lines in Fig. 2, both cases produce at ambient temperature for approximately 20 years.When thermal break-through occurs, the production temperature declines gradually.Over time, the extraction temperature approaches a constant extraction temperature of about 55 • C. The results are identical for both cases because the production rate is scaled proportionally to the reservoir thickness.The solid lines in Fig. 2 show the production temperature for the simulations when heat exchange with the confining layers is enabled.Both cases produce more than 20 years at ambient temperature.When thermal breakthrough occurs, the production temperature decline is less strong due to the effect of heat recharged from the confining layers.The well temperature of the thin reservoir model declines slower compared to the well temperature of the thicker reservoir and also moves towards a higher steady-state extraction temperature.The markers in Fig. 2 represent the production of 1 E rep , and show that as a result of vertical heat recharge more heat is extracted from the geothermal system.This emphasises the positive effect of vertical recharge on heat production from thin reservoirs.This comparison of simulations including and excluding heat exchange between the vertical confining layers and the reservoir highlights the contribution of the heat from the vertical confining layers to the total heat production of a geothermal well.The dotted lines in Fig. 2 correspond to the right hand-side y-axis and indicates the relative difference in production temperature of a geothermal well between the two simulations.These results provide the following insights: (i) The relative effect of vertically recharged heat from the confining layers is larger for a 50-m thick reservoir than for a 200-m thick reservoir, because the ratio of cooled down reservoir area in contact with the confining layers relative to the overall reservoir volume (or production volume) is much larger for this reservoir.
We conclude that generally speaking, for thin reservoirs the relative contribution of the recharge is larger.(ii) The effect of the confining layers on the thermal breakthrough time t BT increases when vertical recharge is taken into account, and is stronger for a 50-m thick reservoir than for a 200-m thick reservoir: it is 21 years for both models excluding the recharge, 23 years for the model with recharge and 200-m thickness, and 28 years for the 50 m-thick model with recharge.The vertically recharged heat from the confining layers also affects the rate at which the well temperature drops, resulting in higher extraction temperatures for thin layer.

The effect of reservoir thickness
Fig. 3 shows the simulation results for varying reservoir thickness with constant production rate and well spacing.The results provide the following insights: (i) To produce the energy content of one representative reservoir volume, the decrease in well temperature in a thinner reservoir is less than in a thicker reservoir.(ii) For thicker reservoirs, the difference between produced heat and vertically recharged heat is larger than for thinner reservoirs.Hence, the relative contribution of vertically recharged heat to the total produced energy is smaller for a relative thicker reservoir.(iii) In thinner reservoirs, the well temperature decreases slower which results in the fact that the representative volume of energy is produced in a 4.2 times longer time period (155/37 years) in a thicker reservoir, while the E rep is only 4 times larger (200/50 m), see Table 2.As a result, the cooled-down reservoir surface area that interacts with the confining layers effectively increases with a smaller reservoir thickness (Fig. 4), in the same period of time.Fig. 4 shows that after 300 years of production, the cooled-down surface area at the top of the reservoir is much larger for a thinner reservoir.The larger cooled down area results in a larger contribution of vertical recharge.After 300 years of production, 81.3 PJ has been vertically recharged for a reservoir of H = 50 m, which is 73% of the energy content, compared to 45.6 PJ, i.e. 58% for H = 200 m.With a relative smaller reservoir thickness the absolute and relative contribution of vertically recharged energy over time is larger.
These results are consistent with the observations in Section 3.1 that with increasing reservoir thickness, the effect of the confining layers becomes less, and as a result, the recharge ratio decreases.During long-  term operation, most of the heat in thin reservoirs is coming from confining layers.

The effect of production rate
Fig. 5 shows the simulation results for a model with H = 200 m and varying production rates.For these simulations, the highest production rate (600 m 3 /h) results in the highest energy output, and thus the fastest cooling of the reservoir per time unit.When comparing the effect of a change in production rate per produced representative volume of energy E rep , the same trend is observed: there is a larger temperature drop per produced E rep when moving towards a higher production rate and a larger difference between the produced and recharged energy.The decrease in recharged heat with an increase in production rate is caused by the speed of the recharged heat.The change in surface area between the representative reservoir volume and the confining layers does not play a role here, because it is equal for all cases at the moment that one representative volume of energy has been produced.It takes 304 years to produce one E rep with Q = 150 m 3 /h in the case of Q = 600 m 3 /h the same amount of energy is produced in 79 years.This means that in the lowest production rate case considered, the energy coming from the vertically confining layers has more than 200 years more time to reach the reservoir.So when taking more time to produce a certain amount of energy with a low production rate, this results in relatively more vertically recharged heat per produced E rep .This effect can be further

Table 2
Schematic overview of the results after production of 1 E rep , the parameters from which the effect is tested and that are varied are indicated in bold.illustrated by analyzing at the vertically recharged heat per year for the cases with Q = 150 and 600 m 3 /h.The yearly vertical recharge increases respectively from 0.26 PJ to 0.99 PJ.Due to the larger cooled down surface area at the interface between reservoir and confining layers, the absolute yearly recharge is larger in the high production rate scenario.This is also caused by the fact that in this short time period the temperature gradient is also large, resulting in a high conduction rate.
For the low production rate scenario the temperature gradient reduces over the longer production time, resulting in a much lower average E rv .Fig. 6 shows the heat vertically recharged relative to the representative reservoir volume ϵ rv against time for the cases with production rate Q = 150 and 600 m 3 /h.It clearly shows that a higher production rate results in a faster vertical recharge.After 50 years of production with a rate of 600 m 3 /h ϵ rv is 0.10.It takes 4 times longer to produce the same amount of energy, with a production rate that is 4 times smaller.Fig. 6 shows that after 200 years of production at a rate of 150 m 3 /h, ϵ rv is 0.19.This means that to produce the same amount of energy a lower production rate allows for more vertical recharge.In the final column of Table 2 the average yearly vertically recharged heat is provided.Because this rate is larger for larger production rates it shows a different trend.This is caused by the fact that at larger Q the E rep is extracted faster (e.g.79 y, compared to 304 y for Q = 600 and 150 m 3 /h, respectively).As a result, at large Q the cooled-down surface area interacting with the confining layers remains at a relatively large temperature difference compared to the temperature in the confining layers.As the heat flux depends on this temperature difference, this results in a larger E rv /y.At low production rates the heat has more time to conduct to the reservoir, resulting in declining temperature difference and thus heat flux.Note that the ratio of the time it takes to produce the E rep with Q = 150 m 3 /h relative to Q = 600 m 3 /h is 3.8 (304 y/79 y), which is smaller than the inverse ratio between the rates, which is 4 (150 m 3 /h/600 m 3 /h), indicating that the absolute larger vertical recharge promotes the yield of the geothermal system, in case of a lower production rate.The results of the varying production rate scenarios provide the following insights: (i) With an increase in production rate, the well temperature declines stronger and the energy content of one representative volume is produced faster.(ii) Despite the much larger absolute recharge at large production rate, the low production which allows the heat to flow from the confining layers into the reservoir yields the largest relative contribution of vertical recharge.

The effect of well spacing
Fig. 7 shows the simulation results for a model with H = 200 m, Q = 150 m 3 /h and varying well spacing.Larger well spacing results in an increase in cooled-down reservoir surface area that interacts with the confining layers, Fig. 8.This stimulates the vertical recharge of heat from the confining layers towards the reservoir, and reduces the rate at which the well temperature drops.The same mechanism discussed in Section 3.2 explains the results of varying well spacing, i.e. large contact surface area with the confining layers promotes vertical recharge.An increase in well spacing results in a longer time before the energy content of the representative volume is produced.In this case, the increased time results in a larger relative and absolute contribution of vertical recharge.

Vertically vs. laterally recharged heat
Due to the cold water injected into the reservoir, heat will flow towards this relatively cool reservoir volume.The heat flow is not only coming from the confining layers, but also laterally, from the reservoir around this cold reservoir volume.Per unit surface area of the interface between cooled down and hot subsurface material, the lateral flow of heat is relatively larger because it is transported both by conduction and free convection, whereas vertical heat exchange is solely a result of conduction, as no convection is possible due to the low permeability in the confining layers.The temperature gradient remains larger at the lateral interface, while at the vertical interface the gradient reduces over time.Despite these two aspects (convection and larger ΔT), the difference between lateral and vertical recharge per unit surface area is small, this is illustrated by the very small difference in location of the markers on the dashed lines in Fig. 2, indicating that the relative contribution of the free convection component of the lateral heat transport is small.While vertically recharged heat will not affect the heat production of neighbouring fields, the laterally recharged heat may interfere with other geothermal sources.For a single geothermal project considered in isolation, the source of the recharged heat is not very relevant, but this is different when they become part of a geothermal network, for example in urban or greenhouse areas.The separation between the lateral and vertical recharge components then becomes important to optimally utilise the geothermal resource and the available heat in its confining layers.
We use the representative reservoir volume (V) to quantify the effect per parameter on the lateral recharge, and examine the amount of energy flowing in both lateral and vertical direction.The sum of these two is the total amount heat inflow, i.e. recharged heat E r .Fig. 9 shows the amount of laterally recharged energy as a function of the total recharged energy content ϵ r , which is determined using Eq.(2.12).A larger inflow of vertically recharged heat then also means a relative lower amount of laterally recharged heat.We find that for the three analysed variables, the lowest amount of laterally recharged heat is obtained with a thin reservoir, a low production rate and a large well spacing.The choice of representative reservoir volume used for these calculations is somewhat arbitrary and a change in shape or size will affect the numeric outcomes.However, the relative differences and the volume are calculated in the same way in all experiments, hence the effect of the exact shape and size Fig. 9. Laterally recharged heat ϵ rl against the total recharged heat ϵ r .The figures show from left to right the effect of reservoir thickness, production rate and well spacing. of this volume cancels out.
In our analysis, we observe effects of the following: (i) reservoir thickness; an increase in reservoir thickness results in an increase of the relative contribution of laterally recharged heat.For a reservoir with H = 200 m, more than 50% of all the energy flowing into the representative volume is coming from laterally confining layers.We relate this to an increase in horizontal surface area of the representative reservoir volume, which allows for more laterally recharged heat.With a well distance of 1200 m and a reservoir of H = 50 m, the ratio of surface area with laterally and vertically confining layers is 1 16 (top and bottom surface area: 2⋅2⋅1200 2 = 5.76 M m 2 , circumference area: 2⋅50⋅1200 + 2⋅2⋅50⋅ 1200 = 0.36 M m 2 , see Fig. 1).While for a reservoir with H = 200 m, this ratio is 1 4 .For well distances of 1000 and 1500 m, these ratios are: 1 3.3 and 1 13.3 and 1 5 and 1 20 , respectively.(ii) production rate; an increase in production rate stimulates the total flow of recharged heat.However, the relative contribution of vertical recharge reduces.This is caused by the fact that, with a higher production rate the temperature difference between the representative reservoir volume and its confining layers decreases faster.This stimulates a higher flow of recharged heat per time unit.Fig. 10 shows the temperature distribution in the reservoir and confining layers after production of one E rep .The figure shows that with a higher production rate, less energy has been produced from the vertical confining layers.Consequently, relative more energy is extracted from the lateral confining layers.This is also a time effect because in the high-productionrate scenario the E rep is produced in much shorter time (Table 2).(iii) well spacing; with a larger well spacing the laterally recharged heat has a smaller relative contribution to the total recharged heat.The increase in surface area with the vertically confining layers allows a higher percentage vertically recharged energy.We should note that the choice of the size of the representative reservoir volume decreases quadratically with the decrease of well spacing.This may cause the cold front to stretch across the borders of the representative reservoir volume, which then contributes to the laterally recharged heat.

Relation to current applications of geothermal energy
The results of this study show that in general the geothermal resource is best utilised when (i) the production rate is low, (ii) the wells are far apart and (iii) the reservoir is thin.This contradicts usual requirements for geothermal projects to maximise production rates.
We observe the following: (i) As discussed in Section 3.3, maintaining a low production rate means that the energy output of the geothermal system per unit time is lower than when a higher production rate is used.If the geothermal project vision is short term, it is obvious that a high production rate is preferable.However, for a more sustainable production, a lower rate will eventually provide a higher energy output, because the confining layers are used more efficiently.(ii) As discussed in in Section 3.4 a larger well spacing increases the energy output of the geothermal project, due to the larger cooleddown surface area that stimulates the recharge of heat.With a well spacing increasing from 1000 to 1500 m, the representative volume (Eq.(2.6)) increases 2.25 times, while the recharged energy (E r ) per produced energy content of a representative volume E rep increases by more than 3 times (from 16 to 48.8 PJ).This shows that a larger well spacing makes more efficient use of the heat available in the confining layers.(iii) Optimising the total heat production from a given reservoir requires a perspective which goes beyond the maximisation of production from a single doublet.To maximise heat production, thereby extracting the maximally possible amount of energy from the reservoir, a low production rate combined with a larger well spacing will be preferred.(iv) A thicker reservoir is preferable due to its higher energy content and larger potential production rate.To cost effectively utilise thin reservoirs, the production such reservoirs will have to be increased to fulfill the current energy demand.The Geological Survey of The Netherlands is currently (2020) examining if thin geothermal reservoirs are economical feasible in the Netherlands.As Section 3.2 illustrates, the contribution of recharged heat will increase the potential of thin reservoirs.Unlike an oil or gas reservoir, in a geothermal project a considerable part of the energy resource is available from the surroundings of the reservoir.Therefore, a shape of the boundary that has a large contact area between the reservoir and its confining layers is preferable as it will stimulate the recharge of heat.
To make more efficient use of geothermal energy it is important to look at (A) larger timescales: beyond the economic lifetime of the geothermal project, and (B) the long term development of the heat demand.Please note from Figs. 4, 7 and 10 that the heat extraction propagates beyond the representative volume.The same representative volume is used by the Dutch Ministry of Economic Affairs and Climate Policy to set the spatial extent of geothermal exploration licences (Mijnlieff and Van Wees, 2009;Mijnlieff, 2014;de Zwart et al., 2012).The license precondition is that the temperature outside of this representative volume boundary does not change within the 30-years' license horizon.This implies that when heat production continues after thermal breakthrough, because production temperature is still useful as a result of heat recharge, these permit limits may need to be reconsidered.In the light of the energy transition, it could be attractive to meet the present energy demand with relatively high rates, resulting in a faster well temperature drop and lower energy output in the future.Technical developments, such as better insulation and more efficient use of the heat, can ultimately lead to a lower temperature demand during the lifetime of a geothermal project.This might increase further efficiency, viability and sustainability of present projects.Better insulation of houses will furthermore lower the temperature demand which could then be supplied by a partially cooled-down reservoir.At the same time, in order to make sustainable use of the subsurface, we could favour the long term use of a geothermal reservoirs and adapt the production strategy accordingly, considering the joint production of multiple doublets rather than optimising per production unit.For further optimisation, it is interesting to explore how multiple, high rate and densely-spaced geothermal projects in a thick reservoir perform, compared to one or two geothermal projects at a low rate in a thin reservoir, utilising the same reservoir volume.Assessing also the financial aspects should provide insight in the financial feasibility of the utilisation strategies presented in this paper.

Whether to include confining layers
In previous research little attention was paid to the contribution of heat coming from the confining layers.It is often neglected or not quantified.Bauer et al. (2019) looked at 2300-m thick reservoir model without considering its environment.Also Crooijmans et al. (2016) did not take any confining layers into account.The fact that the confining layers are often neglected in studies, might be due to time it takes to observe the effect of the confining layers, usually after thermal breakthrough (Fig. 2).And indeed, the performance of the geothermal well before thermal break-through is hardly affected by the confining layers.On the time scale of decades or centuries, the role of the recharge from confining layers becomes more important, as it may extend the lifetime of a geothermal project.Its role can be so large that the total heat recharge may substantially exceed the heat initially in place in the reservoir.Hence recharge from the confining layers should be included for the assessment of total heat recharge.

Limitations/future work
This study focuses on the effect of the confining layers to a geothermal project, using a model that is a simplified representation of a geothermal reservoir as its geology is considered homogeneous and the production takes place with vertical wells.Shetty et al. (2018) and Wang et al. (2020) have demonstrated that the heterogeneity of the reservoir can have a non-negligible effect on the recharge from the confining layers we will consider these aspects in future research.Although our setup captures the main processes, the extent to which our findings are valid in a more realistic conditions remains to be evaluated.The results of our study form a starting point for optimisation strategies when exploring more realistic field conditions.These strategies could consider variations in permeability for placement of inclined or horizontal wells.The cooled-down interface interacting with the confining layers strongly affects the contribution of vertical recharge.An assessment of the placement of the injector to stimulate heat exchange with the confining layers seems a logical extension of this study.
In this study we did not vary any of the subsurface/reservoir material properties.Of course these also affect the rate in which recharge will occur.For generic insights on how recharge can contribute to heat production, further research to identify favourable conditions may be helpful.However, like is the case with the thickness of the reservoir, these conditions cannot be changed.In this study we focus on the parameters that designers have influence on, to optimise the geothermal heat output.In practice, many of these properties are uncertain and spatially variable.It is therefore important for individual projects or field development studies to assess the effect of possible parameter uncertainty and variability on the recharge and heat produced, e.g.Daniilidis et al. (2017).

Conclusion
In this study, we evaluated the vertical heat recharge into a geothermal reservoir from the vertically confining layers.The results show that recharge of heat from the confining layers may contribute considerably to the total energy that a geothermal system produces.This effect is not always taken into account in the design of a geothermal project.The most important effect of the contribution of heat recharge is that the production temperature after thermal breakthrough may exhibit a considerably smaller decline, resulting in a longer lifetime of the geothermal project.The largest contribution of heat recharge to the heat production can be expected for systems with: • thin layers: The relative effect of vertically recharged heat from the confining layers is largest for a thinner reservoir, because the ratio of cooled-down reservoir area in contact with the confining layers relative to the overall reservoir volume (or production volume) is much larger for a thinner reservoir.Hence, for thin reservoirs the relative contribution of the recharge is larger.Of course the total amount of heat produced from a thin reservoir is smaller than from a thick reservoir.• low production rate: With an increase in production rate, the well temperature declines stronger and the energy content of one representative volume is produced sooner.Despite the much larger absolute recharge at large production rate, the low production rate allows the heat to flow from the confining layers into the reservoir and yields the largest relative contribution of vertical recharge.• large well spacing: Larger well spacing results in an increase in cooled-down reservoir surface area that interacts with the confining layers.This stimulates the vertical recharge of heat from the confining layers towards the reservoir, and reduces the rate at which the well temperature drops.
In general, the results indicate that a more efficient use of the available geothermal heat requires operational conditions which are currently not always met.This is due to the relative short time horizon used to assess (financial) feasibility of geothermal projects, while longer time horizons are needed to be able to take the benefit from heat recharge into account.As the thermal recharge that result from these conditions increase the lifetime and the total heat production of a geothermal project considerably, they will ultimately lead to a better use of available geothermal resources.Building on the present work, the optimisation of well placement and reservoir-wide production and licensing strategies helps to promote vertical heat recharge and to optimise field utilisation for multiple geothermal systems.This work may serve as a stepping stone for a more comprehensive evaluation of the economic lifetime of geothermal projects and the long-term viability of geothermal energy supply.
.1) with ρ 0 the fluid density at reference concentration and temperature [kg/m 3 ], ρ s is the density of the source and sink [kg/m 3 ].μ 0 is the dynamic viscosity at reference concentration and temperature [kg/m/s].K 0 is the hydraulic conductivity tensor of material saturated with reference fluid [m/s].h 0 is the hydraulic head measured in terms of the reference fluid of a specified concentration and temperature [m].S s,0 is the specific storage, defined as the volume of water released from storage per unit volume per unit decline of h 0 [m 1 ].t is time, C is salt concentration [kg/m 3 ] and q ′ .7) in which c b is the volumetric bulk heat capacity of the reservoir [J/m 3 / K], V is the representative volume [m 3 ], T av is the average initial reservoir temperature [ • C] and T inj is the injection temperature [ • C].
Fig. 1.A schematic representation of the model (thickness H [m]) illustrates the reservoir with confining layers.Injection and production well are separated by well distance (L).The representative volume (V) in the middle of the reservoir has dimensions H⋅L⋅2L.

Fig. 2 .
Fig. 2. The effect of confining layers, and its vertically recharged heat on the well temperature.On the left y-axis the well temperature for two reservoir (H = 50 & 200 m) with and without the confining layers during 300 years of production.On the right y-axis is the deviation between the well temperature with and without the confining layers modeled.The markers represents 1 E rep that has been produced.
E.A.M. de Bruijn et al.

Fig. 3 .
Fig. 3.The effect of reservoir thickness on recharge.The produced and vertically recharged energy content against the production temperature.Measured for a reservoir of H = 50, 100 & 200 m during 300 years of production with a constant production rate of 300 m 3 /h for all cases.Each marker on the line represents 100 years of production.

Fig. 4 .Fig. 5 .
Fig. 4. The temperature distribution at the bottom of the top confining layer after 300 years of production, with on the left a reservoir of 50 m thick, on the right a reservoir of 200 m thick.The black contour line indicates the 89 • C boundary.The red triangle pointing down and up represent the injection and production well, respectively.The orange box indicates the representative volume boundary.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 6 .
Fig. 6.The ϵ r,v during a production period of 300 years for the cases with Q = 150 and 600 m 3 /h.

Fig. 7 .
Fig. 7. Production temperature against a dimensionless axis.Measured for a reservoir of H 200 m during 300 years of production, with a production rate of 150 m 3 / h, and a well spacing of L = 1000, 1200 and 1500 m.Each marker on the lines represents 100 years of production.

Fig. 8 .
Fig. 8. Temperature distribution in the layer overlying the reservoir, after production of 1 E rep , with on the left side the model with L = 1000 m, and on the right side the model with L = 1500 m.The black contour line indicates the 89 • C boundary.The red triangle pointing down and up represent the injection and production well, respectively.The orange box indicates the representative volume boundary.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 10 .
Fig. 10.Temperature distribution in a vertical cross section of the reservoir and confining layers after production of one E rep .With on the right the model with Q = 150 m 3 /h and on the left the model with Q = 600 m 3 /h.The red lines represent the wells, and the black lines the upper-and lower bound of the reservoir.The orange box gives the boundaries of the representative reservoir volume boundary.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Table 1
Model input parameters.