European Geologist Journal 43

Assessment of deep geothermal energy potential in Northern and Central Portugal

By M.M. Miranda*, N.V. Rodrigues, Willis-Richards and A.J.S.C. Pereira

* CEMMPRE – Centre of Mechanical Engineering Materials and Processes, Department of Earth Sciences, Rua Sílvio Lima Polo II, University of Coimbra. 3030-790 Coimbra, Portugal,


Northern and Central Portugal are being intensively studied to identify geological environments with high potential for Enhanced Geothermal Systems/Hot Dry Rock (EGS/HDR). The Beira granite is the selected target for geothermal feasibility studies. This granitic mass shows a heat production of 5 µW.m-3. Surface heat flow is estimated to be 100–120 mW.m-2. The depth for a possible geothermal reservoir suitable for electricity generation is predicted to be 6 km. Simulations carried out to assess the potential for the development of an EGS reservoir on the Beira granite show that, although these reservoirs in Portugal are not as easily created by stimulation as they might be in other regions, reservoirs of potential commercial interest were generated in the majority of model realisations.

1. Introduction and background geology

Northern and Central Portugal are characterised by wide outcrops of Hercynian granitic rocks Middle to Upper Carboniferous in age emplaced into a metasedimentary sequence, the Schist-Greywacke Complex, of Pre-Cambrian to Lower Palaeozoic age. The granitic rock masses were emplaced during and after the last thickening and subsequent post-collisional extensional events of the Hercynian orogeny (D3 ductile compressive deformation and D4 brittle extensional phases) that occurred from the Mid-Carboniferous to Permian. Their emplacement was associated with a high temperature event that occurred mainly during and after the ductile deformation phase D3 (Namurian–Westphalian in age).

These regions of mainland Portugal may have significant potential for enhanced geothermal system style geothermal energy. This is supported by the occurrence of several thermal springs with temperatures ranging up to 76 ºC (IGM, 1998) in the outcropping Hercynian granites. The available heat flow density maps (e.g., Correia and Ramalho, 2005) suggest that these regions have surface heat flow values in the range 65–80 mW.m-2, with an average of 78 mW.m-2. However, there are no reliable deep temperature measurements from the granitic plutons themselves to support this conclusion.

In the present work we discuss the geological evidence that suggests higher heat flow values are likely and that pinpoints, within the several granitic units, which one has the highest geothermal potential for Enhanced Geothermal Systems/Hot Dry Rock (EGS/HDR) exploitation. Studies of the heat-producing elements (HPEs) and their host minerals, rock thermal properties and estimates of radiogenic heat production (RHP) were carried out, leading to estimates of temperature-at-depth and simulations of an EGS reservoir within the target area.

2. Heat-producing elements, radiogenic heat production and rock thermal properties

A total of 314 samples were collected from the different granitic units outcropping in northern and central regions (Figure 1). The methods used to assess the HPEs concentrations and the results obtained are described in detail in the work of Lamas et al. (2017). Briefly, it was observed that there is a substantial increase of U and Th concentrations from the oldest granitic units to the younger ones, as can be seen in Figure 2. The A-group (A1 and A2) shows a U average concentration of 4.5 and a Th average concentration of 9.1 (Lamas et al., 2017). For the C-group the U average concentration is 5.5 for C1, 8.2 for C2, and 7.5 for C3 (Lamas et al., 2017). Regarding Th average concentrations in the C-group, C1 shows a value of 22.9, C2 21.3, and C3 shows an average Th concentration of 18.7 et al., 2017).

Figure 1: Geographical and geological framework with location of the samples analysed and the different granitic units studied. A to C indicate each granitic unit, with A being the oldest granites emplaced in the northern and central regions, and C the younger ones; “Samples” refer to the surface samples, and “AQ1 borehole” to the subsurface samples. For more explanation see Lamas et al. (2017).

Figure 2: Tukey boxplots showing the U and Th distribution amongst the different granitic units studied. Adapted from Lamas et al. (2017).

As U and Th concentrations tend to increase from the older granitic units to the younger ones, the same tendency is observed for radiogenic heat production, as discussed in Miranda et al. (2015b) and shown in Figure 3. The A-group presents a surface heat production of 2 µW.m-3 while C-group 4 µW.m-3 (Miranda et al., 2015).

Figure 3: Tukey boxplots showing the surface radiogenic heat production distribution amongst the different studied granitic units. Adapted from Miranda et al. (2015).

Based on the highest U and Th concentrations, and consequently the highest radiogenic heat production, and due to their estimated thickness (5 km mean and 10 km maximum; Machadinho, 2014) the C2-lithology was chosen as the target lithology for more detailed studies. Lamas et al. (2015a) pointed out that within the C2-unit, the granitic rock mass that occupies the central region of Portugal (between Viseu, Guarda and Almeida; Figure 1), known as the Beira granite, is the most radiothermal granite.

The existence of a deep well near Almeida (see Figure 1) with core recovery allowed us to study the distribution of the heat-producing elements with depth (see Lamas et al., 2017, for further details). Due to hot spring fluid flow from fractures, this well is not suitable for estimation of heat flow. The U concentration at Almeida increases from surface to subsurface by about 30 %. At the surface its average concentration is 9.9 (Lamas et al., 2015), whilst in the subsurface it increases to 14.3 (Lamas et al., 2017). Consequently, the heat production of the Beira granite increases with depth. The surface value is 4.5 µW.m-3 while at depth it increases to 5.2 µW.m-3 (Miranda et al., 2015).

From our study of the subsurface samples it was shown that the radioelements, and consequently the heat production, have a remarkably uniform distribution along 1,000 meters depth (Figure 4; Lamas et al., 2017), a trend noted in other deep boreholes of the Hercynian chain (e.g., Sams and Thomas-Betts, 1988 and references therein).

Figure 4: Heat-producing elements (U, Th and K2O) and radiogenic heat production (A) distribution at depth.

Our study of the subsurface samples highlights that even when surface samples are apparently freshest, the mobility of U due to weathering had taken place, and thus surface samples are not representative of the whole rock mass for geothermal studies.

Beyond the heat-producing elements and heat production distribution, we also studied the thermal properties of samples from the Beira granite. The method used and the results obtained are described in Miranda (2014). The average thermal conductivity obtained is 3.1W.m-1.K-1, for thermal diffusivity the obtained value is 1.6 x10-6 m2.s-1, and thermal capacity 1.9 x106 J.m-3.K-1 (equivalent to 720 J.Kg-1.K-1). Their distribution with depth is shown in Figure 5.


Figure 5: Distribution of thermal properties with depth, Almeida borehole. λ – thermal conductivity; k – thermal diffusivity, ρcp – thermal capacity.

3. Heat flow and temperature modelling

Heat flow and geothermal gradient estimations were carried out for the Beira granite using an approach described in Pereira et al. (2016). A 100 km linear traverse was selected within the granite of interest, in such a way that AQ1 borehole at Almeida (see Figure 1 and Figure 6) was the central point of the model.


Figure 6: Detail of the geological map seen in Figure 1 with the location of the linear segment studied for heat flow and temperature estimations. From Pereira et al. (2016).

A conceptual crustal model for the area of interest was built based on available geophysical and geological data (see Pereira et al., 2016 for a detailed discussion), and was assumed to be as depicted in Figure 7. The heat production was assumed to follow a stepwise distribution in which each crustal layer is composed by a uniform value. In this way, the values of 5.2, 10.9, 1.4, and 0.8 µW.m-3 were considered for the granite, upper crust II, mid-crust, and lower crust, respectively. The value of heat production assumed for the low-velocity layer upper crust II is based on a purely mathematical hypothesis and it was calculated by the vp-A relationship, in which A increases with decreasing vp(e.g., Rybach and Bunterbath, 1984).

Figure 7: Conceptual crustal model for the studied segment. From Pereira et al. (2016).

A similar assumption of each crustal layer being composed by uniform values was followed for thermal conductivity, in which the granitic layer was considered to have a value of 3.1W.m-1.K-1, the upper crust II and mid-crust the value of 2.5 W.m-1.K-1, and the lower crust was given the value of 2.1 W.m-1.K-1. For the schist layer thermal conductivity of 3.0 W.m-1.K-1 and heat production of 1.6 µW.m-3 were considered (Pereira et al., 2016).

The simulations of present-day in-depth temperatures and heat flux distributions were carried out using MATLAB, based on the numerical solution of the 2D steady-state heat conduction equation over the rectangular domain shown in Figure 7. A uniform Moho heat flux (q”) of 25 mW.m-2, characteristic in Palaeozoic folded units (Cermak, 1993), was assumed to be entering the domain through the lower horizontal boundary and the top boundary is considered to be at a constant and uniform surface temperature of 16 ºC. The vertical boundaries are both assumed to be adiabatic (i.e., no horizontal heat flow) (see Pereira et al., 2016 for more details).

Different model geometries and several thicknesses for a schistose metasedimentary cover were assessed, enabling us to study the effect of each layer thickness and the effect of a thermal blanket on heat flow and temperature. For the granitic part of the model the surface heat flow was estimated to vary from 114 to 130 mW.m-2, and the geothermal gradient from 32.5 to 37.3 ° (Pereira et al., 2016). The average surface heat flow obtained is 123 mW.m-2, which is higher than previously published synoptic data (e.g., Correia and Ramalho, 2005). For the geothermal gradient the average value is 35.3 ° These values are good approximations; however, they are dependent on the distribution model used for heat production and thermal conductivity. Taking into account such uncertainties, but assuming reasonable confidence in extrapolating measured subsurface heat production to depth in the Beira granite, we can conclude that the surface heat flow in the Beira granite is higher than has been proposed to date.

The modelled heat flow and temperature curves, showing the effect of the metasedimentary cover, are shown in Figure 8. Both profiles were modelled only to 6 km depth since this is the likely maximum acceptable depth for a potentially commercially viable EGS project.


Figure 8: Temperature and heat flux curves. From Pereira et al. (2016).

Assuming a metasedimentary layer with 4 km thickness the surface heat flow is estimated to be around 100 mW.m-2; however, with a thinner cover, the surface heat flow reaches an estimated value of near 115 mW.m-2 (Pereira et al., 2016). At 6 km depth, the average temperature in the granitic part of the modelled segment is 210 ºC, decreasing to 205 ºC assuming a metasedimentary cover of 1 km and to 185 ºC if the metasedimentary cover has a thickness of 4 km (Pereira et al., 2016). These temperatures are within the range of temperatures for potential EGS reservoirs.

4. Numerical modelling of EGS reservoir development in Beira granite

As discussed above, due to temperature values found, the Beira granite is the most promising scenario for an EGS reservoir in Central Portugal. Numerical simulations to assess its potential for the development of such reservoirs were carried out using the FRACSIM3D model, as is discussed in Miranda et al. (2017).

FRACSIM3D is a 3-D stochastic fracture network and rock mechanics model that addresses the geology and geomechanics that control the development of an EGS reservoir created by hydraulic stimulation of crystalline hot rocks. A detailed description of this model is present in the work of Jing et al. (2000) and the version used here is updated in respect to fracture dilation laws and the introduction of a cohesion term to mimic the strengths of compressional asperities that have to break before slip takes place. This allows crude estimates of potential microseismic event magnitudes. Briefly, FRACSIM3D uses approximations to follow the complex changes that occur within the rock mass during stimulation, and it simulates the subsequent steady state circulation of the heat exchange system created, thus helping in the design of a reservoir for the situation at hand. The model is written as a sequence of modules which, starting from available field measurements, simulates a fracture network in a cube-shaped model volume. It follows the stimulation of fractures from one and more wells at a given injection pressure and then calculates water circulation through the shear-dilated, compliant fracture network. Flow rates, water losses, heat extraction, tracer particle tracking and a frictional law allowing both slip strengthening or weakening are implemented. In summary, FRACSIM3D can analyse both the creation/ enhancement of the geothermal reservoir as well as the approximate operation of it once it has been developed.

As inputs, FRACSIM3D requires data regarding:

  • rock elastic and thermal properties;
  • tectonic stresses;
  • fracture orientations, frequency, and fracture length distribution;
  • fracture friction, roughness, closure curve, macro-asperity strength expressed as a cohesion term;
  • location and properties of any possible major faults;
  • undisturbed large scale permeability;
  • in-situ fluid pressure (hydrostatic);
  • rock temperature at depth.

Based on these inputs and well orientations and positions, stimulation pressures, rock volume to be treated, and the injection and recovery well pressures during circulation, it is possible to obtain as model outputs:

  • injection and recovery flow rates;
  • thermal performance over time;
  • tracer response curve;
  • maps of flow, pressure, transmissivity, and temperature;
  • maps of fracture slip and unstable slip;
  • synthetic well flow logs.

Our imperfect knowledge of the input parameters (e.g., in-situ large scale permeabilities can only be guessed, and might be set at levels 10x greater than for other EGS systems, thus representing the worst case scenario for water loss) leads to high uncertainties in the modelling process. Nevertheless, not giving too much emphasis to the absolute values of the flow rates, with progressive improvements in reservoir stimulation design and circulation conditions, our Beira granite scenario shows some quasi-commercial flow rate EGS model reservoirs (Table 1), bearing in mind that a commercial EGS reservoir is expected to have flow rates higher than 10 l/sec. With the treatment options studied (reduction of circulation pressure and increase of stimulation pressure in both wells, with addition of a small extra stimulation of the recovery well) water losses were reduced from 29 l/sec to 18 l/sec. Although this still is a relatively high water loss rate, stimulation design and circulation treatments must be assessed.

A comparison study was also carried out between Beira granite and a Rosemanowes (Cornwall) scenario at the same depth using the same in-situ permeability, conservative stimulation pressures of 1 MPa below the minimum effective in situ principal stress, and assuming asperities with strength sufficient to generate micro-seismic events. The differences between the two models were the stress field orientation, fracture orientations, and fracture length distributions. Results obtained from a number of realisations are shown in Table 2.

Table 1: Flow rates obtained by changing circulation and stimulation pressures. Inter-well distance of 300 m; open-hole lengths of 600 m; -30º as the rotation applied to well system; -40º as the rotation applied about the y-axis and 20º as the rotation applied about the z-axis for each well mid-point. The azimuth of the maximum horizontal principal stress was considered 140º.

Table 2: Flow rates for Rosemanowes and Beira granitic units.

The Beira scenario shows a mean recovery flow of 95 %, while for Rosemanowes it was 99 %, with flows 54 % greater. This performance difference is due mainly to the difference in fracture orientations relative to the principal stresses. To generate equivalent flows in Beira, the stimulation pressure has to be increased, and either a second smaller high pressure stimulation of the recovery well or proppants must be introduced, or possibly both.

From multiple realisations and sensitivity studies carried out using the Beira granite scenario (to assess uncertainty and potential investment risk), we conclude that EGS reservoirs in granites in Portugal are not as easily created by stimulation as they might be in the Rhine Graben (France, Germany) or in SW England. This is mostly due to the rotation of the Iberian Peninsula in Cretaceous times that has resulted in few joints being optimally oriented for slip at the present day. These fractures, not very well-oriented for fluid flow, reduce flow by about 50 %. Necessary stimulation pressures are estimated to be about 4–5 MPa higher than those that would be required at similar depths in Soultz-sous-Fôrets or Cornwall. Nevertheless, we generated reservoirs of potential commercial interest in the majority of model realisations. Further work is required to investigate the existence of possible deep faults or hydrothermally altered zones as potential EGS reservoirs; such features may be common and may have in-situ apertures, frictional and shear dilation properties distinct from the joint-based scenario reported above.

5. Conclusions

The granitic units cropping out in Northern and Central Portugal were studied to assess their heat production potential. Amongst these, the Beira granite showed highly radiothermal heat production, and was thus selected as a target for geothermal investigations regarding EGS reservoirs.

Although we have very imperfect knowledge of the input parameters for the numerical modelling of an EGS in the Beira granite, these preliminary results show that EGS reservoirs are not as easily created by stimulation as they might be in other regions. However, reservoirs of potential commercial interest were generated in the majority of the simulations carried out.


This work has been framed under the Initiative Energy for Sustainability of the University of Coimbra and supported by the Energy and Mobility for Sustainable Regions Project (CENTRO/07-0224-FEDER-002004), co-funded by the European Regional Development Fund (ERDF) through the «Programa Operacional Regional do Centro 2007 – 2013 (PORC)», in the framework of the «Sistema de Apoio a Entidades do Sistema Científico e Tecnológico Nacional», and by the «Fundação para a Ciência e Tecnologia». The authors would like to acknowledge Prof. José António Simões Cortez and Câmara Municipal de Almeida for allowing the study of the borehole core.


Cermak, V. 1993. Lithospheric thermal regimes in Europe. Physics of the Earth and Planetary Interiors, 79. 179–193.

Correia, A., Ramalho, E.C. 2005. Updated surface heat flow density map in Mainland Portugal. Proceedings, World Geothermal Congress 2005. 1–5.

IGM 1998. Recursos Geotérmicos em Portugal Continental: Baixa Entalpia (Geothermal resources in Mainland Portugal: Low-enthalpy systems). Available online at:

Jing, Z., Willis-Richards, J., Watanabe, K., Hashida, T. 2000. A three-dimensional stochastic rock mechanics model of engineered geothermal systems in fractured crystalline rock. Journal of Geophysical Research, 105. 23663–23679.

Lamas, R., Miranda, M.M., Pereira, A.J.S.C., Ferreira, N., Neves, L.J.P.F. 2015. Distribuição dos elementos radiogénicos nas rochas granitóides aflorantes na Zona Centro-Ibérica (Centro e Norte de Portugal) (Distribution of radiogenic elements of the granitic rocks outcropping in the Central Iberian Zone (Central and Northern Portugal)). Proceedings, X Congresso Ibérico de Geoquímica. 173–176.

Lamas, R., Miranda, M.M., Pereira, A.J.S.C., Neves, L.J.P.F., Ferreira, N., Rodrigues, N.V. 2017. 3-D distribution of the radioelements in the granitic rocks of Northern and Central Portugal and geothermal implications. Journal of Iberian Geology, DOI: 10.1007/s41513-017-0001-y.

Machadinho, A. 2014. Modelação da geometria de rochas granitóides recorrendo a métodos geofísicos gravimétricos e magnéticos: uma contribuição para a avaliação do potencial geotérmico na região Centro de Portugal (Modeling of the geometry of granitic rocks using gravimetric and magnetic geophysical methods: a contribution for the assessment of the geothermal potential of Central Portugal). PhD dissertation, University of Coimbra..

Miranda, M.M. 2014. Propriedades termofísicas em rochas graníticas: o caso da sondagem profunda de Almeida (Guarda, Portugal central) (Thermophysical properties of granitic rocks: the case of the deep borehole of Almeida (Guarda, central Portugal)). Master thesis, University of Coimbra.

Miranda, M.M., Lamas, R., Pereira, A.J.S.C., Ferreira, N., Neves, L.J.P.F. 2015. Potencial térmico das rochas granitóides aflorantes na Zona Centro-Ibérica (Portugal) (Heat production of the granitic rocks outcropping in the Central Iberian Zone (Portugal)). Comunicações Geológicas, 102(Especial 1). 133–136.

Miranda, M.M., Willis-Richards, J., Rodrigues, N.V. 2017. Numerical modelling of engineered geothermal systems (EGS) using FRACSIM3D in granitic rocks of central Portugal. Proceedings, 8th European Geothermal PhD Day.

Pereira, A.J.S.C., Costa, J.J., Panão, M., Miranda, M.M., Machadinho, A., Lamas, R., Neves, L.J.P.F., Rodrigues, N.V. 2016. Estimation of heat flow and geothermal gradient from numerical modelling in central Portugal. Proceedings, European Geothermal Congress 2016. 8 pp.

Rybach, L., Bunterbath, G. 1984. The variation of heat generation, density and seismic velocity with rock type in the continental lithosphere. Tectonophysics, 103. 335–344.

Sams, M.S., Thomas-Betts, A. 1988. 3-D numerical modelling of the conductive heat flow of SW England. Geophysical Journal, 92. 323–334.

This article has been published in European Geologist Journal 43 – Geothermal – The Energy of the Future