New methods of geothermal potential assessment in the Pannonian basin

Abstract

The Pannonian basin in Central Europe is well known for its rich geothermal resources. Although geothermal energy has been utilised, mainly for direct use purposes, for a long time, there are still a lot of untapped resources. This paper presents novel methods for outlining and assessing the theoretical and technical potential of partly still unknown geothermal reservoirs, based on a case study from the Dráva basin, one of the sub-basins of the Pannonian basin along the Hungarian–Croatian border. The presented methods include reservoir delineation based on combining geological bounding surfaces of the Upper Pannonian basin-fill units with a set of isotherms deriving from a conductive geothermal model. The geothermal potential of each identified reservoir was calculated by a Monte Carlo method, which was considered as being represented by the heat content of the fluids stored in the effective pore space (‘moveable fluid’). The results underline the great untapped geothermal potential of the Dráva basin, especially that of the reservoir storing thermal water of 50–75°C, which has the largest volume and the greatest stored heat content.

Introduction

The Pannonian basin in Central Europe is one of the European areas with a well-known positive geothermal anomaly, and whose rich geothermal resources have long been utilised mainly for direct use purposes. Hungary was ranked as the fifth country in Europe in 2017 for geothermal district heating and thermal-water heating cascade systems, with estimated amounts of 223.36 MW_t installed capacity and 635.66 GW_h production (EGEC, 2018). Individual space heating (mostly associated with spas) adds a further estimated installed capacity of about 77.2 MW_t and 83.1 GW_h annual production (2017), whilst the key player is still the agriculture sector, especially in southern Hungary, where heating of greenhouses and plastic tents accounts for about 358 MW_t installed capacity (Nádor et al., 2019).

Despite non-technical barriers (complex and time-consuming licensing processes, deficiency in capital investment), the geothermal sector in Hungary is developing fast, with an average of 10–20 new geothermal wells drilled each year. This adds to the outstandingly high number (about 980) of already existing thermal water wells in Hungary (defined in national legislation as wells with outflow temperature higher than 30°C), with a production history of often 30+ years. The number of emerging projects demands a firm geoscientific basis for assessment of the available geothermal resources, not only at local (project-scale) level, but also at regional and even at trans-boundary scales, as co-use of the spatially extended geothermal aquifers might have serious impacts on thermal groundwater flow systems (Rman et al., 2016; Szőcs et al., 2018) carrying the subsurface heat medium, especially in regions where production is extensive (e.g. southern part of the Great Hungarian Plain).

Project developers need to have realistic numbers on the available geothermal resources. Unfortunately, as pointed out by Rybach (2013), the generic term ‘geothermal potential’ is commonly used without a clear qualification on the specific type of potential (Table 1), which might confuse investors as to what performance can be expected from a given geothermal prospect.

Table 1. Categories of geothermal potential (Rybach, 2010)

Although thermal karst systems in the deep-lying Mesozoic carbonates in the basement store a significant amount of thermal water in Hungary (Korim, 1972; Goldscheider et al., 2010), the majority of thermal water production comes from thick Neogene basin-fill clastic sequences. Therefore in the present paper we focus on the geothermal potential assessment methods targeting these porous reservoirs of the Pannonian basin. We present delineation methods for the so-called ‘basin-fill’ (BF) reservoirs storing thermal water of different temperatures which have already been applied to different parts of the Pannonian basin (Rotár-Szalkai et al., 2017), and also discuss novel methods for quantifying the extractable amount of heat associated with the heat content of fluids stored in the effective pore space (i.e. ‘moveable water’). We present a case study from the Dráva basin extending from SW Hungary to N Croatia, where a probabilistic approach was applied to quantify the exploitable amount of heat from various potential porous reservoirs. This case study is part of a larger project – DARLINGe http://www.interreg-danube.eu/approved-projects/darlinge – which targets the sustainable and more efficient use of the existing yet still largely untapped geothermal resources at the southern part of the Pannonian basin.

Geology and geothermal conditions of the Pannonian basin

The heat flow density of the Pannonian basin ranges from 50 to 130 mW m⁻² with a mean value of 90–100 mW m⁻² (Hurter & Haenel, 2002; Lenkey et al., 2002; Horváth et al., 2015). The wide range of values is explained by the presence of recharge areas and the cooling effect of the infiltrating meteoric waters. The overall positive geothermal character of the Pannonian basin with a geothermal gradient of about 45°C km⁻¹ (Dövényi and Horváth, 1988) is related to the Early–Middle Miocene crustal extension. As a result of mid-Miocene rifting (Horváth and Royden 1981; Horváth 1993), deep basins originated with the deposition of thick offshore pelagic sediments in the rapidly sinking basin interiors, and coarse-grained brackish littoral facies sediments at the margins and on the tilting basement highs. This was followed by a significant thermal subsidence since the Late Miocene (Royden et al., 1983), when a huge inland lake (Lake Pannon) with a brackish water body formed about 11 Ma as a remnant of the Middle Miocene Paratethys Sea (Magyar et al., 1999; Harzhauser and Piller, 2007). This giant lake basin was gradually filled up by forward accretion of sediment packages during the Late Miocene and Early Pliocene, originating in the surrounding uplifting Alpine and Carpathian mountain belts. The resulting siliciclastic sequences with low thermal conductivity (the ‘Pannonian’ sedimentary succession) can be as thick as several thousand metres in the deepest sub-basins (Figures 1 and 2). The ‘Lower Pannonian’ sedimentary package consists of coarse material on basin floors deposited by turbidity currents with subsequent sedimentation of silt and argillaceous marl on the slopes. Prograding shelf sedimentation resulted in the deposition of thick sand bodies in the shelf front environment, overlain by fine-grained sediments of the shelf and alluvial plain, representing the latest stages of the filling-up of the lake basin, referred to as ‘Upper Pannonian’ (Bérczi and Phillips, 1985; Magyar et al., 1999, 2013). The most favourable porous geothermal aquifers are formed in these ‘Upper Pannonian’ shelf front sands at depths from 500 to 2000 m, where the temperature ranges from 40 to 110°C (Szanyi and Kovács 2010; Nádor et al., 2012; Rman, 2014; Horváth et al., 2015).

Figure 1. Thickness of Neogene sediments in the Pannonian basin (Horváth, 1988). DB refers to the Dráva basin case study area, the top and bottom surfaces of whose most important potential porous geothermal reservoirs are displayed in Figures 3–6, also showing the location of the regional geological cross section in Figure 2.

Figure 2. Geological cross section of the Dráva basin (Rotár-Szalkai et al., 2018). Location is shown in Figure 1.

Under the several kilometres thick Neogene sedimentary sequences the basement of the Pannonian basin is extremely diverse and shows a complex structure, built up of various metamorphic and non-metamorphic Palaeo- and Mesozoic crystalline and carbonate rocks whose formation was associated with the Alpine–Carpathian orogene. They are arranged into nappes along thrust sheets, dissected by strike-slip and normal faults, associated with the multiphase tectonic development of the basin (Fodor et al., 1999; Csontos and Vörös 2004; Horváth et al., 2006; Schmid et al., 2008; Haas et al., 2012; Haas & Budai, 2014). These rocks also form some outcropping mountain regions, so-called ‘inselbergs’. The permeability of some basement carbonate rocks is enhanced by karstification that affected them during their geological evolution. The crystalline basement rocks might also serve as potential geothermal reservoirs, especially along regional tectonic zones with increased fissure permeability, or in their uppermost weathered-altered zones of 15–30 m average thickness. At the basement depth (around 2000 m or more) temperature can exceed 100–150°C (Dövényi et al., 2002).

Methods of reservoir delineation

The method for delineating ‘potential reservoirs’ within the several thousand metre thick Pannonian (especially Upper Pannonian) sedimentary succession has been elaborated previously (Rotár-Szalkai et al., 2017). It is based on the combination of representative geological horizons (top and bottom bounding surfaces of a given hydrogeological unit, potentially storing thermal water) by various subsurface isotherm maps by applying SURFER and ArcGIS softwares. As the regional thermal water flow systems of the Pannonian basin often cut across national borders, editing of such geological and geothermal model horizons requires massive data harmonisation among countries sharing the same transboundary geothermal reservoirs (Maros et al 2012; Nádor et al., 2012; Tóth et al., 2016; Szőcs et al., 2018).

Despite the lower permeability of the clayey–marly strata within the thick shelf front to alluvial plain sedimentary sequences of the Pannonian basin-fill complex, hydraulic connection exists between the regionally extended sandy aquifer layers with a hydraulic conductivity of about 4 × 10⁻⁶ to 5 × 10⁻⁵ m s⁻¹. This makes the entire sedimentary series one large hydrostratigraphic unit (i.e. a unit encompassing different rock formations but possessing similar hydrogeological features) characterised by an almost uniform hydrostatic pressure.

The further subdivision of this large hydrostratigraphic unit into potential basin-fill geothermal reservoirs was done in three consecutive steps. First the bottom surface of the shelf-front sandy deposits and the top surface of the fine-grained alluvial sediments as bounding geological horizons of the main geothermal aquifers were edited.

Than various isotherm surfaces were selected to represent different temperature ranges of various utilisation schemes important in this region: thermal waters with temperatures of 30–50°C are mainly used for balneology; temperatures between 50 and 100°C are primarily suitable for direct heat utilisation, whereas the 75–100°C subcategory represents the temperature range in which thermal water can be economically used for larger district heating systems. Reservoirs having temperatures above 100°C have the potential for combined heat and power generation (CHP) projects using binary technologies. The various subsurface isotherm surfaces were edited on the basis of a simplified conductive model (details discussed in Rotár-Szalkai et al., 2018), which was established for the entire southern part of the Pannonian basin in the DARLINGe project. The basic hypothesis of the calculation method was that the higher heat flow (higher values than world average) in the Pannonian basin is caused by the thinning of the lithosphere that resulted from lithosphere stretching coeval with the formation of the basin (Royden et al., 1983; Dövényi & Horváth 1988). The relief of the pre-Pannonian basement evolved parallel with the intensive inflow of sediments due to thermal subsidence, i.e. the depth of the basin is proportional to the degree of thinning of the crust. The spatial variation of heat-flow density therefore reflects changes in the basin depth. Calculations of the depth value of the isotherm surfaces (grids) in the basin-fill sediments were performed by multiplying the basin depth by the average geothermal gradient. To validate the modelled temperature distribution, we compared them with the real temperature measurement data from wells. The comparison showed very close statistical agreement, but locally some differences occurred. Differences are expected in the regions where convection (modifying effect of groundwater flow) has an important role in spatial subsurface distribution of temperature, or is caused by uncertainty of temperature measurements.

As a last step the different basin-fill reservoirs were spatially outlined by combining the grids of the geological bounding horizons and the respective isotherms. For example, the combination of the geological horizon representing the top of the the fine-grained alluvial sediments and the 50°C isotherm surface designates the minimum depth of available thermal water with temperature higher than 50°C (top of BF50–75 reservoir). Similarly, the combination of the same geological horizon and the 75°C isotherm surface defines the minimum depth of available thermal water with temperature higher than 75°C (top of BF75–100 reservoir) and so on. It has to be highlighted that coming from the nature of sedimentary units, the basin-fill reservoirs are lens-shaped 3D bodies, whose top and bottom bounding horizons are bent surfaces which join each other, therefore they do not have ‘side walls’. Therefore the top (plan) view of a reservoir visualises the spatial distribution of both the top and bottom bounding surfaces.

Methods for potential assessment

When assessing the technical potential of a given geothermal reservoir, we want to express in a quantitative way the amount of geothermal energy to be potentially exploited. The first limit is the maximum realisable drilling depth which is 5–7 km (Muffler & Cataldi, 1978). It is also important to provide reliable values for the ‘recovery factor’ (the ratio of heat recoverable versus heat in place) that are needed to convert theoretical potentials into technical potentials. According to practical experience, its value is rather uncertain and varies in wide ranges. Depending on the geology and permeability, it is between 0.1 and 0.25 in porous systems and between 0.08 and 0.2 in fractured systems in the Pannonian basin (Nádor, 2016).

Since in large sedimentary basins, such as the Pannonian basin, the extraction of heat from the rock generally occurs by the production of thermal water (carrier of the heat energy), only the heat content of the ‘moveable fluid’ (fluid situated in the effective pore space) provides realistic values of the exploitable heat amount. In our study we considered this heat content as equivalent to the theoretical potential (physically usable energy as defined by Rybach, 2010) and performed calculations based on the classical formula of heat-in-place (i.e. volumetric method) (Muffler and Cataldi 1978), such as the following:

(1) $${Q_{{\rm{effpor}}}} = C \times m \times {\rm \Delta} T$$

• C = 4,187,106 J/(1000 kg·°C) is the specific heat capacity of water

• m = mass of moveable water (effective porosity) in reservoirs (kg)

• ΔT = T _average – T _reference is the excess temperature over 30°C (°C).

In our study we paid special attention to the determination of the effective porosity (which stores the moveable water that can be exploited), therefore the mass of moveable water in the reservoirs was calculated according to:

(2) $$m = \rho \times {\rm\Phi_{{\rm{eff}}}} \times {\rm{ }}{\rm \Delta} H \times A$$

• ρ = 1.0 [1000kg m⁻³] = specific density of water

• Φ_eff = average effective porosity of the reservoir (V_{eff pore}/V_total, no unit)

• ΔH = average thickness of the reservoir (m)

• A = area of the reservoir (m²)

where the effective porosity is:

(3) $${{\rm \Phi} _{{\rm{eff}}}} = (1 - {\rm{SH}}){\rm{ }} \times {{\rm \Phi} _{\rm{t}}}$$

• SH = average clay content of the reservoir (V _clay/V _total, no unit)

• Φ_t = average total porosity of the reservoir (V _pore/V _total, no unit).

Effective porosity is derived from total porosity. Its value is the function of actual basin depth (effect of compaction) and of the clay content. The clay content in our model was handled as laminar clay (the lithology consists of 100% sand and 100% clay layers, no disseminated clay). Lacking detailed geological information on the exact facies distribution of the Upper Pannonian shelf front, shelf and alluvial plain sequences in the southern part of the Pannonian basin, the clay content in the model was considered as varying between 20 and 60% (based on expert estimation and literature of the previous regional geological and hydrogeological modelling (Tóth et al., 2016), with an average value of 40% for the entire sequence.

Results: Dráva basin case study

The Pannonian basin, as a large sedimentary basin complex, is composed of several sub-basins (Figure 1). The Dráva basin, situated along the border of Hungary and Croatia (Figure 1), is such a sub-basin, where the Neogene sedimentary succession is up to 5–7 km in total thickness (Figure 2) and is composed of three second-order megacycles separated by four major erosional unconformities (Stafic et al., 2003).

In this very thick basin-fill complex we focused only on the Upper Pannonian – Pontian shelf front, shelf-plain, coastal plain to alluvial plain sedimentary succession (Pannonian Lake nearshore facies (BF) in Figure 2), as this is the major thermal-water bearing unit of the Dráva basin, which accumulated in an aggradational style due to the persistent high rate of subsidence coupled with an enormous rate of sediment supply. The top and bottom bounding horizons of this sedimentary unit – as characteristic facies, thus lithological boundaries – could be mapped from well logs and seismic sections, which were partly available for the study area (Horváth et al., 2012). The various isotherms (subsurface temperature maps) of 30, 50, 75, 100, 125 and 150°C derived from the conductive geothermal model. The top and bottom surfaces of the different basin-fill reservoirs were than edited by combining the geological bounding horizons and the respective isotherms. The top and bottom surfaces of the most important reservoirs of the Dráva basin (BF50–75, i.e. reservoirs storing thermal water with temperature higher than 50°C, and BF75–100, i.e. reservoirs storing thermal water with temperature higher than 75°C) are shown in Figures 3, 4, 5 and 6.

Figure 3. Top surface of the BF50–75 reservoir. The white patch indicates that the 50°C isotherm is most probably found in the very thick Quaternary sediments in the central depression of the Dráva basin (Figure 2), which is not considered as a sensu stricto ‘basin-fill reservoir’ due to its relatively low temperature.

Figure 4. Bottom surface of the BF50–75 reservoir.

Figure 5. Top surface of the BF75–100 reservoir.

Figure 6. Bottom surface of the BF75–100 reservoir.

Once the reservoirs were spatially delineated, the next step was to calculate the effective porosity heat content and the recoverable heat for each reservoir category according to Equations 1, 2 and 3, summarised in Table 2. The distinctive parameters of a geothermal reservoir (area, thickness, porosity, and reservoir temperature) vary within a certain range of values, which can best be assessed by probabilistic approaches. Quantification of uncertainties in the parameters of the probability distribution can be dealt with quite well using the Monte Carlo simulation method (e.g. Shaha et al., 2018). The Monte Carlo simulation is a stochastic method, which uses randomness to solve problems that might be deterministic in principle. The computerised mathematical technique relies on repeated random sampling of input data series.

Table 2. Input and calculated parameters for quantifying the technical potential

In the Dráva basin case we applied the Monte Carlo simulation during the calculations of the effective porosity (Equation 3), where an average value and a standard deviation for the error of the estimated total porosity (±10%) and for the estimated clay content (±20%) were applied. For the rest of the parameters, static values were used: the area of the reservoirs was determined as the surface projection of the top of the reservoirs. The thicknesses of the reservoirs vary from place to place, and were individually determined in a 500 × 500 m grid from the geological model horizons as the distance between the top and bottom horizons. Temperature was estimated from the conductive model as an average for the given reservoir considering the depth and the thickness. These input parameters are summarised in Table 3.

Table 3. Input parameter values for the Dráva basin

Heat energy stored in the effective pore space (Equation 1) was also calculated by the Monte Carlo method where all the possible results were obtained with a cumulative distribution function in the form of P90, P50 and P10 cases representing high, medium and low levels of confidence. Two types of estimation were done: calculating the total heat content of moveable fluids (theoretical potential), and calculating the technical potential considering the recovery factor with a conservative value of 0.1 (Table 4). An example of the distribution of theoretical potential according to the various probability distributions is shown for one selected reservoir (BF75–100°C) in Figure 7, while Figure 8 summarises the calculated total heat content of moveable fluids for all of the BF reservoir types in the Dráva basin.

Table 4. Estimated heat content of effective porosity, at confidence levels P90, P50, P10, with and without applying the recovery factor of 0,1

Figure 7. Distribution of the estimated theoretical potential according to P90, P50 and P10 for the BF75–100 reservoir in the Dráva basin.

Figure 8. Calculated heat content of the BF reservoirs in the Dráva basin. Note that the BF50–75 reservoir is the most propitious in terms of P10, P50 and P90 values.

Discussion and conclusions

The paper has presented novel approaches to geothermal resource estimation at regional scales in two aspects: (1) it introduced methods for delineating potential geothermal reservoirs at basin scales, even extending beyond the political borders of neighbouring countries, and (2) it highlighted the importance of the extractable amount of geothermal energy stored in the effective porosity, and presented a probabilistic approach for its quantification.

(1) In the case of large, extended transboundary geothermal aquifers with multiple users having cumulative impacts by their thermal water abstractions, it is of utmost importance to have a joint and regional understanding of the reservoir dimensions and properties irrespective of state borders in order to guarantee sustainable production levels. The Pannonian basin in Central Europe is a typical case, where clarifications of the same reservoirs vary in the neighbouring countries due to the use of different geological nomenclatures for the same formation, different understandings of the thermal water flow systems, dissimilar models of the subsurface temperature field, etc. Although some promising progress in terms of common models has been achieved in the western part of the Pannonian basin (e.g; Rman et al. 2016; Rotár-Szalkai et al., 2017; Nádor et al., 2019), the southern part remained a great challenge due to the highly uneven data distribution, and lack of FAIR principles (Findable, Accessible, Interoperable and Reusable) for data policies (Nádor, 2019). Nevertheless through the joint efforts of geoscientists from the DARLINGe project countries (Hungary, Slovenia, Croatia, Serbia, Bosnia and Herzegovina, and Romania) it was possible to edit some harmonised geological and geothermal model horizons (the latter based on a simplified conductive model), which served as the basic input information to delineate common, transboundary geothermal reservoirs within the several thousand meters thick Neogene sedimentary succession of the Pannonian basin. The simple but efficient method of reservoir delineation was based on the spatial combination of the bounding geological surfaces of the main geothermal aquifers with the respective isotherms, which made it possible to define a series of reservoir categories i.e. porous sediment bodies that store thermal water in temperature ranges of 30–50°C (BF 30–50), 50–75°C (BF 50–75), 75–100°C (BF 75–100), etc.

(2) The second group of new results is associated with the resource assessment of the delineated reservoirs. Estimations of geothermal resources – even for areas with well-acknowledged geothermal potential, such as the Pannonian basin – are often exaggerated. To provide scientifically based and reliable information, it has to be emphasised that potential investors should consider only the extractable amount of geothermal energy, which depends on many technical and non-technical limiting factors. In the case of sedimentary basins, one such aspect is the effective pore space, which stores the moveable (thermal) groundwater, i.e. the carrying medium of geothermal energy to be exploited. While the total porosity depends on a series of controlling factors, such as the original sediment composition and structure, grain size, etc., all determined by the former depositional environment, the effective porosity is rather governed by diagenetic processes, and its value is the function of actual basin depth (effect of compaction) and of the clay content.

All these novel methodological approaches were tested and presented through a case study of the Dráva basin, part of the Pannonian basin along the Hungarian–Croatian border. Figure 8 and Tables 3 and 4 let us draw the conclusion that both in volume and stored heat content, the BF50–75 reservoir (i.e. the part of the sandy Upper Pannonian sedimentary succession that stores exploitable thermal water with temperature between 50 and 75°C) is the most important in the Dráva basin. Its top surface is found around −800 to −900 m a.s.l. (Figure 3), representing an easy and relatively cheap drilling target. Although the BF30–50 reservoir is almost as extensive in volume (1,669,707 km³) as the BF50–75reservoir (1,819,788 km³), it has only about half as much stored heat content, due to its shallower depth (Table 4). The BF75–100 reservoir, situated at an average of −1400 to −1500 m depth (Figure 5), is more restricted in dimensions (associated with the central depression of the Dráva basin), but its heat content is still significant (about half of the BF50–75 reservoir’s; Table 4). These numbers underline the region’s great potential for space and district heating, which at the moment exists only at a few locations. Even though reservoirs storing thermal water above 100°C are limited to relatively small areas and bigger depths (below 2000 m), they represent a potential for CHP applications.