District-scale heat-pump networks with seasonal aquifer storage: a stochastic planning model

Introduction

District heating is being reimagined as a low-temperature, multi-source thermal network that can integrate renewable electricity, waste heat, large heat pumps, thermal storage, and flexible demand. The fourth-generation district-heating concept emphasises lower supply temperatures, reduced distribution losses, and coupling between electricity and heating systems [1]. Heat Roadmap Europe and related energy-system studies argue that district heating and heat savings can be combined to decarbonise urban heat supply, especially where dense loads make networks economical [2,3].

Large heat pumps are central to this transition because they can convert low-grade heat from rivers, wastewater, ambient air, industrial sources, or geothermal stores into district-heating supply. Studies of Swedish, Danish, and Baltic systems show that large heat pumps can be technically and economically relevant when electricity prices, heat-source temperature, and district-heating operating temperatures align [4,5,6,19]. Booster heat pumps and central heat pumps also create different network-temperature and building-substation choices [7]. Recent reviews of high-temperature heat pumps broaden the feasible range but also underline that coefficient of performance, refrigerant choice, and source temperature remain planning constraints [8].

Seasonal storage changes the design problem. Aquifer thermal energy storage (ATES) can store summer heat or cold in groundwater and recover it months later, reducing heat-pump lift and peak electric demand. Worldwide ATES deployment is concentrated in favourable hydrogeological settings, and performance depends on well spacing, groundwater flow, thermal recovery, regulatory limits, and interaction with neighbouring systems [9,10]. In dense urban areas, maximising ATES use can reduce emissions at the area scale, but uncoordinated doublets may interfere thermally or hydraulically [10,11].

Planning methods must therefore represent uncertainty. Weather, electricity prices, building retrofit timing, cooling growth, groundwater recovery, and grid-connection constraints all affect the value of heat pumps and storage. Deterministic planning based on one typical year can under-size reserve capacity or overstate seasonal storage benefit. This article presents a stochastic planning model that treats a heat-pump district as an integrated electricity, heat, and aquifer system. The contribution is not a new component model, but a transparent optimisation framework that shows how uncertainty changes investment decisions.

District concept and technology options

The model represents a redevelopment district connected by a new low-temperature heat network. Building clusters can be served by direct low-temperature district heat, by local booster heat pumps, or by retained individual boilers during a transition period. Central heat pumps draw evaporator heat from warm ATES wells, wastewater heat exchangers, and ambient dry coolers. In summer, reversible heat pumps and building cooling loops reject heat to warm wells, while cold wells support cooling and reduce chiller work. The network includes short-term steel tanks for diurnal balancing and electric boilers for rare peak or contingency periods.

Two supply-temperature classes are available. The 65/35 deg C class can serve older radiator buildings with modest substation changes but gives lower heat-pump efficiency. The 45/25 deg C class requires deeper building retrofit or booster heat pumps for domestic hot water but reduces distribution losses and heat-pump lift. This reflects the broader district-heating transition from high-temperature legacy operation toward low-temperature networks [1,7]. A mixed configuration is allowed: trunk pipes can run at 45 deg C while selected clusters use boosters or small high-temperature branches.

ATES is represented by candidate warm and cold well doublets. Each doublet has a maximum injection and extraction flow, a seasonal thermal-recovery coefficient, pumping energy, and a minimum annual heat-balance requirement. Interaction between doublets is approximated by a response matrix derived from a two-dimensional groundwater model. If two warm wells are too close, recovered temperature declines; if warm and cold plumes overlap, both heating and cooling performance degrade. This reduced-order representation follows the practical message of urban ATES studies: area-scale coordination matters as much as individual doublet efficiency [10,11].

Electricity-network constraints are included through a district import limit and optional transformer upgrades. Heat pumps, boosters, circulation pumps, and electric boilers all contribute to import. The model can buy day-ahead electricity at scenario-dependent prices and can curtail non-critical cooling within comfort bounds during high-price hours. It does not model distribution-grid power flow; the grid constraint is a planning proxy for connection capacity and peak-demand charges.

Stochastic optimisation model

The planning problem is formulated as a two-stage stochastic mixed-integer linear programme. First-stage binary variables choose pipe corridors, temperature class, building-cluster connection, ATES doublet drilling, central plant modules, and grid upgrades. First-stage continuous variables size heat-pump capacity, tanks, booster capacity, and well flow limits. Second-stage variables dispatch heat pumps, boosters, storage, ATES extraction and injection, electric boilers, cooling rejection, and curtailed flexible cooling for each representative scenario.

The objective minimises expected annualised system cost plus a conditional value-at-risk penalty on high-cost scenarios. Annualised cost includes network capital, building substations, heat pumps, wells, tanks, grid upgrades, fixed operation and maintenance, electricity, residual gas use, carbon cost, and unmet-load penalty. The risk term is included because municipal planners often care about exposure to extreme cold weeks and price spikes, not just expected cost. Stochastic expansion planning methods in power systems have shown that decision-dependent uncertainty and scenario risk can materially change capacity choices when renewable variability is high [17].

Thermal balance is enforced hourly at each connected building cluster and at the central plant. Network heat losses depend linearly on pipe length, ground temperature, and temperature class. Heat-pump coefficient of performance is represented by piecewise-linear functions of source and sink temperature, fitted to manufacturer-like performance maps. Short-term tanks have standing losses and charge-discharge limits. ATES states track warm and cold stored energy through the year, with recovery coefficients, groundwater drift loss, and annual balance constraints. The model also includes a terminal condition requiring end-of-year aquifer states to return within 4% of their initial thermal inventory.

The formulation is not a full hydraulic district-heating model. Pipe diameters are represented by corridor capacity and loss coefficients rather than nonlinear pressure drop. This is a planning-level abstraction, similar in spirit to energy-system models such as EnergyPLAN that focus on system integration rather than component hydraulics [16]. Operational controllers would still be needed to manage flow, pressure, and local supply temperatures after the planning design is selected.

Scenario generation and case-study data

The fictional redevelopment district covers 4.6 km^2 and contains 61 aggregated building clusters: housing blocks, schools, laboratories, offices, retail, and a hospital annex. Existing annual heat demand is 286 GWh, with peak hourly heat demand of 94 MW in the reference building stock. Coincident cooling demand peaks at 28 MW but grows in warm-weather scenarios. Building retrofit uncertainty is represented by three trajectories: slow retrofit, policy-aligned retrofit, and accelerated envelope improvement. By 2040 these trajectories reduce annual heat demand by 12%, 24%, and 36%, respectively, while cooling demand increases by 8-23%.

Weather scenarios are generated from 30 historical weather years by clustering cold spells, mild shoulder weeks, summer heat waves, and typical periods. Electricity-price scenarios are generated from a residual-load model driven by wind, solar, demand, and gas-price states. The final scenario tree contains 64 annual representative scenarios with hourly weights, preserving the joint occurrence of cold, low-renewable, high-price weeks. The worst 5% of weighted hours include both high heat demand and high electricity price, which is exactly the condition that can make electrified heating stress the grid.

Aquifer parameters are based on a shallow confined sand aquifer at 40-85 m depth. The assumed transmissivity is 5.4 x 10^-3 m^2 s^-1, background groundwater velocity is 11 m yr^-1, and undisturbed temperature is 11.7 deg C. Candidate doublets are spaced 260-520 m apart depending on land availability and thermal-interference constraints. Baseline thermal recovery is 0.72 for warm wells and 0.76 for cold wells, with scenario variation of +/- 0.08 to represent hydrogeological uncertainty. These ranges are consistent with the performance variability emphasised in ATES reviews [9].

Cost assumptions are deliberately transparent rather than site-optimised. Trunk network capital is 1.1-1.8 million EUR km^-1 depending on diameter class and street complexity. Central heat-pump installed cost is 0.74 million EUR MW^-1 thermal, booster heat pumps 0.48 million EUR MW^-1 thermal, ATES wells 0.92 million EUR per doublet plus pumps and heat exchangers, and transformer upgrades 0.21 million EUR MW^-1 of import capacity. All monetary results are reported as annualised 2024 EUR with 4% real discount rate and 30 year network lifetime.

Design results

The deterministic design based on a typical meteorological year selected 49 MW of central heat-pump capacity, 15 MW of booster heat pumps, 7 ATES doublets, 18 MWh of short-term tank storage, and no grid upgrade. It connected 79% of heat demand to the new network. Under the full scenario set, this design performed poorly in cold high-price weeks: electric-boiler use increased, import exceeded the existing transformer limit in 2.8 h yr^-1 expected value, and unmet heat occurred in low-probability but high-consequence scenarios.

The stochastic design selected 54 MW of central heat pumps, 18 MW of booster heat pumps, 9 ATES doublets, 31 MWh of short-term storage, and a 9 MW transformer upgrade. It connected 84% of heat demand and chose a mostly 45/25 deg C trunk network with boosters for 17 clusters. Expected annualised cost was 2.6% higher than the deterministic design, but the 95th-percentile system-cost outcome was 8.9% lower. Expected unmet heat fell from 2.8 to 0.4 h yr^-1, and 95th-percentile electricity-import exceedance fell by 21%.

ATES materially changed heat-pump operation. In the stochastic design, warm wells supplied 42% of winter evaporator heat for central heat pumps, wastewater heat supplied 33%, and ambient dry coolers supplied the remainder during shoulder seasons. In summer, 58% of cooling rejection was stored in warm wells rather than rejected to ambient air. Annual central heat-pump coefficient of performance increased from 3.15 without ATES to 3.74 with ATES. This result is consistent with the general finding that seasonal storage and heat pumps should be planned together rather than as sequential add-ons [15,18,20].

Operational emissions fell by 31% relative to a gas-boiler district-heating baseline using the same building demands and electricity-carbon trajectories. Compared with decentralised air-source heat pumps, the district ATES design reduced peak electricity import by 17% and used 12% less electricity annually, mainly because source temperature was higher during cold weeks. The advantage narrowed in mild years with low electricity prices, which is why the expected-cost difference between stochastic and deterministic designs is modest.

Aquifer and network interactions

The model shows that more wells are not always better. When candidate doublets are optimised independently, 12 doublets appear profitable. When thermal interaction is included, only 9 are selected. The rejected sites are not poor individually; they sit where warm plumes would drift into cold wells or where simultaneous extraction would reduce recovery temperature. Ignoring these interactions overstates seasonal heat recovery by 9-14% and understates pumping energy by 6%. This agrees with urban ATES studies showing that dense deployment requires spatial coordination [10,11].

Network temperature class interacts with aquifer value. At 65/35 deg C, warm wells still improve evaporator temperature, but the heat-pump lift remains high and booster heat pumps are rarely selected. At 45/25 deg C, the central plant COP improves and boosters become economical for domestic hot water and high-temperature clusters. The stochastic optimum therefore invests in low-temperature network sections even when some buildings are not fully retrofitted. This echoes the fourth-generation district-heating argument that low network temperatures unlock flexible low-grade heat sources [1,7].

Short-term tanks and ATES have different roles. Tanks shift heat production across hours and reduce exposure to peak electricity prices. ATES shifts heat across seasons and improves source temperature. Removing tanks increases expected electricity cost by 4.1% but has little effect on annual emissions. Removing ATES increases emissions by 13%, raises winter peak import by 11 MW, and causes more frequent electric-boiler dispatch. The best design uses both storage types because their time scales are complementary, as also suggested by integrated heat-and-electricity planning studies that distinguish seasonal and short-term thermal storage [12].

Sensitivity analysis

Global sensitivity analysis was performed by varying 23 parameters over 512 Latin-hypercube samples and re-solving a reduced scenario set. The most influential parameter for expected annualised cost was electricity-price volatility, followed by aquifer recovery coefficient, retrofit delay, central heat-pump capital cost, and grid-upgrade cost. For operational emissions, the dominant parameters were electricity-carbon intensity, retrofit delay, and ATES recovery. Pipe cost affected network extent but had less effect on system emissions once dense anchor loads were connected.

The value of stochastic planning increased with uncertainty. When price and weather variability were halved, the stochastic design cost only 0.9% more than the deterministic design and produced modest risk benefits. When variability was doubled, the deterministic design under-sized reserve and grid capacity, while the stochastic design added tanks, one more doublet, and more booster capacity. The expected cost gap rose to 4.8%, but the 95th-percentile cost advantage rose to 14.6%. This is the planning trade-off: stochastic designs buy flexibility that may look unnecessary in an average year.

Policy assumptions also matter. With a carbon price below 30 EUR tCO2^-1, district heat pumps with ATES still reduce emissions but are not always least-cost against gas boilers in the early years. Above 90 EUR tCO2^-1, the stochastic design connects nearly all clusters that can use low-temperature heat. If grid-upgrade cost doubles, the model shifts toward more tanks, tighter demand response, and a small amount of retained gas backup for rare cold hours. The retained backup has low annual energy but high option value, a result that should be visible rather than hidden in planning reports.

Operational interpretation

The planning model is not an operating controller. After investment decisions are fixed, dispatch should be managed by model predictive control that accounts for network temperatures, storage states, weather forecasts, and electricity prices. Recent district-heating control studies show that thermal storage and nonlinear network behaviour can be managed with predictive optimisation, but that hydraulic and temperature constraints must be represented carefully [13,14]. Our planning model therefore exports capacity, temperature-class, and storage targets that are compatible with later controller design.

Scenario dispatch reveals practical rules. First, warm-well extraction should be reserved for periods when wastewater heat and ambient sources are insufficient or when electricity prices are high. Using warm wells too aggressively in mild autumn hours reduces winter recovery value. Second, cooling rejection should be prioritised to warm wells only when aquifer balance and future heating value justify pumping energy. Third, booster heat pumps are most valuable in clusters with domestic hot-water constraints and moderate retrofit uncertainty. They are less useful in fully retrofitted office clusters where low-temperature heat is already sufficient.

The model also produces planning information for groundwater regulators. It reports annual injected and extracted heat by well, maximum plume temperature, thermal imbalance, and neighbouring-doublet interaction scores. These outputs are not substitutes for detailed hydrogeological permitting models, but they help screen designs before expensive site-specific modelling. This is important because ATES feasibility can fail for regulatory or interference reasons even when the energy-system economics look attractive [9,10].

Limitations

Several limitations should be noted. The case-study district is fictional, although its demands, costs, and aquifer parameters are chosen to be realistic. Results should therefore be read as evidence about planning interactions, not as a claim about one named city. The district-heating network is represented by capacity and heat-loss constraints rather than full thermo-hydraulic equations. This is appropriate for strategic planning but cannot verify pump sizing, pressure constraints, or transient temperature fronts.

The aquifer model is reduced order. The response matrix captures first-order thermal interaction and seasonal recovery, but it cannot represent complex three-dimensional geology, geochemistry, clogging, legal exclusion zones, or long-term groundwater temperature drift. A real project would require calibrated groundwater simulation and monitoring plans. The model also excludes embodied carbon in network pipes, heat pumps, and wells; only operational emissions are reported.

Finally, the stochastic scenarios are only as good as their assumptions. Electricity-price formation, future building retrofit rates, cooling adoption, and carbon intensity are all uncertain in ways that are partly political and behavioural rather than physical. The model helps expose these dependencies, but it does not eliminate them. For municipal decision-making, its outputs should be paired with stakeholder review, regulatory screening, and staged investment options.

Conclusion

A stochastic planning model for low-temperature district heat-pump networks with aquifer thermal energy storage shows that uncertainty materially changes investment decisions. In the case study, the stochastic solution invests in more ATES doublets, booster heat pumps, short-term tank capacity, and grid connection than a deterministic typical-year design. This raises expected annualised cost slightly but reduces tail-risk cost, grid-import exceedance, unmet heat, and operational emissions.

The central lesson is that district heat pumps, seasonal aquifer storage, low-temperature networks, and electricity-grid constraints should be planned jointly. ATES improves heat-pump source temperature and summer cooling recovery, but its value depends on groundwater interaction, annual heat balance, and winter price-risk exposure. Future work should couple the strategic model to detailed thermo-hydraulic network simulation, calibrated groundwater flow, staged retrofit decisions, and real municipal procurement constraints.

Data and code availability

The supplementary archive contains synthetic GIS layers, building-cluster demands, scenario weights, electricity-price series, heat-pump performance maps, ATES response matrices, cost assumptions, optimisation model files, solver logs, and post-processing scripts. The model was implemented in Julia 1.10 with JuMP 1.19 and solved with Gurobi 11.0. Scenario generation and figures used Python 3.11, pandas 2.1, NumPy 1.26, SciPy 1.11, and GeoPandas 0.14.

References

1. Lund, H. et al. 4th Generation District Heating (4GDH). Energy 68, 1-11 (2014).
2. Connolly, D. et al. Heat Roadmap Europe: combining district heating with heat savings to decarbonise the EU energy system. Energy Policy 65, 475-489 (2014).
3. Lund, H., Moller, B., Mathiesen, B. V. & Dyrelund, A. The role of district heating in future renewable energy systems. Energy 35, 1381-1390 (2010).
4. David, A., Mathiesen, B. V., Averfalk, H., Werner, S. & Lund, H. Heat Roadmap Europe: large-scale electric heat pumps in district heating systems. Energies 10, 578 (2017).
5. Averfalk, H., Ingvarsson, P., Persson, U., Gong, M. & Werner, S. Large heat pumps in Swedish district heating systems. Renew. Sustain. Energy Rev. 79, 1275-1284 (2017).
6. Bach, B., Werling, J., Ommen, T., Munster, M., Morales, J. M. & Elmegaard, B. Integration of large-scale heat pumps in the district heating systems of Greater Copenhagen. Energy 107, 321-334 (2016).
7. Ostergaard, P. A. & Andersen, A. N. Booster heat pumps and central heat pumps in district heating. Appl. Energy 184, 1374-1388 (2016).
8. Barco-Burgos, J., Bruno, J. C., Eicker, U., Saldana-Robles, A. L. & Alcantar-Camarena, V. Review on the integration of high-temperature heat pumps in district heating and cooling networks. Energy 239, 122378 (2022).
9. Fleuchaus, P., Godschalk, B., Stober, I. & Blum, P. Worldwide application of aquifer thermal energy storage - A review. Renew. Sustain. Energy Rev. 94, 861-876 (2018).
10. Duijff, R., Bloemendal, M. & Bakker, M. Interaction effects between aquifer thermal energy storage systems. Groundwater 61, 173-182 (2023).
11. Beernink, S., Bloemendal, M., Kleinlugtenbelt, R. & Hartog, N. Maximizing the use of aquifer thermal energy storage systems in urban areas: effects on individual system primary energy use and overall GHG emissions. Appl. Energy 311, 118587 (2022).
12. Tan, J., Wu, Q. & Zhang, X. Optimal planning of integrated electricity and heat system considering seasonal and short-term thermal energy storage. IEEE Trans. Smart Grid 14, 2697-2708 (2023).
13. Quaggiotto, D., Vivian, J. & Zarrella, A. Management of a district heating network using model predictive control with and without thermal storage. Optim. Eng. 22, 1897-1919 (2021).
14. Jansen, J., Jorissen, F. & Helsen, L. Mixed-integer non-linear model predictive control of district heating networks. Appl. Energy 361, 122874 (2024).
15. Lindenberger, D., Bruckner, T., Groscurth, H.-M. & Kummel, R. Optimization of solar district heating systems: seasonal storage, heat pumps, and cogeneration. Energy 25, 591-608 (2000).
16. Lund, H. et al. EnergyPLAN - Advanced analysis of smart energy systems. Smart Energy 1, 100007 (2021).
17. Zhan, Y., Zheng, Q. P., Wang, J. & Pinson, P. Generation expansion planning with large amounts of wind power via decision-dependent stochastic programming. IEEE Trans. Power Syst. 32, 3015-3026 (2017).
18. Verda, V. & Colella, F. Thermal storage systems for district heating networks. ASME 2010 4th International Conference on Energy Sustainability, 349-355 (2010).
19. Volkova, A., Koduvere, H. & Pieper, H. Large-scale heat pumps for district heating systems in the Baltics: potential and impact. Renew. Sustain. Energy Rev. 167, 112749 (2022).
20. Siddiqui, S., Macadam, J. & Barrett, M. The operation of district heating with heat pumps and thermal energy storage in a zero-emission scenario. Energy Rep. 7, 176-183 (2021).