Geophysically guided well siting at the Aluto-Langano geothermal reservoir

Ethiopia is rich in geothermal resources, which formed due to active magmatism along the Main Ethiopian Rift (MER). The MER is the northern part of the larger East African Rift, which is a present-day prime example of magma-assisted continental rifting. In the MER, rifting-related magmatic processes lead to magma transport from the lower crust to the surface, where the melt is stored at upper crustal depths of approximately 5 km (Corti, 2009; Samrock et al., 2021). These shallow intrusions are ideal geologic settings for the formation of magma-driven geothermal systems, which are attractive candidates for baseload renewable electricity generation (Scott et al., 2015; Jolie et al., 2021). We note that, in 2010, geothermal power plants installed in the top 10 of the most productive volcano-hosted geothermal fields produced $5536MWe$⁠, which is more than half of the globally installed geothermal power generation capacity of $10,621MWe$ at the time (numerical data from Bertani, 2012). This fact illustrates the remarkable potential of high-temperature geothermal fields to provide energy for electricity generation. Expanding the utilization of the abundant geothermal energy resources in countries along the East African rift, where geothermal energy usage is still in the early stages of development, could provide a significant share of baseload electric power. The socioeconomic impact could be vast, especially in countries such as Ethiopia, where only 44% of the population has continuous access to electricity (Burnside et al., 2021; Gebremeskel et al., 2021). Increasing geothermal development also could mitigate conflicts over Nile ownership in connection with hydro-power projects and downstream water needs (Sterl et al., 2021; Tayia et al., 2021).

Of all countries along the East African Rift, geothermal energy generation from rifting-related magma-driven geothermal resources is presently the most advanced in Kenya. In 2020, Kenya reached an overall installed geothermal power generation capacity of $1193MWe$⁠, of which $599MWe$ have been added since 2015, illustrating an upsurge in geothermal energy utilization with ambitious plans to reach $5000MWe$ by 2030. In Ethiopia, geothermal energy utilization is under development but at a less advanced stage. Presently, 22 high-temperature geothermal systems have been identified within the scope of a geothermal development program that began four decades ago (Kebede, 2016; Hochstein et al., 2017). The identified prospects can all be classified as magma-driven geothermal systems suitable for electricity generation with a total estimated power generation capacity of more than $10,000MWe$ (Kebede, 2016; Burnside et al., 2021). However, only approximately five of the prospects are developed or currently under development.

Geothermal energy utilization in Ethiopia is most advanced at the Aluto-Langano geothermal prospect, which also hosts the country’s first producing geothermal power plant. Aluto is a silicic volcanic complex located between lakes Ziway and Langano (Figure 1) which underwent several eruptive cycles with melt recharge into shallow, upper crustal magma reservoirs (Hutchison et al., 2016). Such shallow magma reservoirs drive powerful hydrothermal convection systems that potentially form prospective geothermal reservoirs (Jolie et al., 2021). Several resemblant volcanic complexes that host geothermal systems are found in magmatic segments along the rift-aligned extensional Wonji fault belt that accommodates 5 mm/year of extension between the Nubian and Somalian plates (Corti, 2009; Kebede, 2016). Other volcanic geothermal prospects in the vicinity of Aluto are Tulu Moye and Corbetti, which are located approximately 50 km north and south of Aluto (Gíslason et al., 2015; Samrock et al., 2018), respectively.

The first productive well at Aluto was drilled in 1983 and reached 2144 m depth, where it encountered maximum temperatures of 320°C (Hochstein et al., 2017). It took until 1999 for a $7.3MWe$ pilot power plant to be installed at the site, which produced energy mainly from two wells, LA-3 and LA-6, drilled in 1983 and 1984, respectively (Hochstein et al., 2017). Although these two early wells have been productive without ceasing since they were drilled, the power plant installed subsequently produced electric energy only with interruptions and mostly at a limited capacity of about half of the nominal capacity of $7.3MWe$ due to other, mostly technical, limitations (Hochstein et al., 2017; Burnside et al., 2021). Currently, expansion work and drilling of new wells have been expected to reach $75MWe$ production in two phases (Cherkose and Mizunaga, 2018; Burnside et al., 2021). Two new wells, drilled in 2013, were productive and could generate 2.7–$3.8MWe$ each (Sisay, 2016; Burnside et al., 2021), whereas three new wells completed in 2021 and 2022 are being tested (see Table 1). The location of the wells is shown on the map in Figure 1. Productive wells are all drilled along a north-northeast-striking fault that intersects Aluto called the Artu-Jawe fault zone (AJFZ). The AJFZ was identified as the major upflow zone of the geothermal reservoir (e.g., Hutchison et al., 2016; Jolie et al., 2019). Another, less prominent, fault is the Worbota-Adonsha fault zone (WAFZ) crosscutting through Aluto northeast of the AJFZ with a strike angle of northeast (Mulugeta et al., 2021) (see also Figure 1). The WAFZ has not been considered for drilling until now, but as we will show, it holds potential for productive drilling.

Numerous geoscientific studies at the Aluto focused on geothermal exploration and the underlying volcanic system in general. Studies embrace a large set of disciplines including geologic mapping (e.g., Hutchison et al., 2016), geochemical analyses and petrological studies (e.g., Gianelli and Teklemariam, 1993; Teklemariam et al., 1996; Pürschel et al., 2013; Gleeson et al., 2017; Jolie et al., 2019), GPS measurements and satellite remote sensing (e.g., Birhanu et al., 2018; Albino and Biggs, 2021), and geophysical methods such as passive seismics (e.g., Nowacki et al., 2018; Wilks et al., 2020), gravity (e.g., Saibi et al., 2012), magnetic mapping (e.g., Mulugeta et al., 2021), and electrical resistivity sounding (e.g., Samrock et al., 2015; Hochstein et al., 2017; Cherkose and Mizunaga, 2018; Samrock et al., 2021). A detailed review of the studies conducted at Aluto is beyond the scope of this study. Here, we focus on the magnetotelluric (MT) imaging of the Aluto volcano, which has been demonstrated to have a high value in geothermal exploration by providing key information on subsurface structures of volcanic magmatic-hydrothermal systems at depths relevant for exploration and academic studies (Younger, 2014; Trainor-Guitton et al., 2017; Trainor-Guitton, 2020).

The MT method is a natural-source geophysical imaging method with virtually no environmental footprint. MT is used to map the 3D electrical conductivity distribution in the subsurface. In geothermal exploration, electrical conductivity is particularly useful because it serves as a proxy for key geothermal parameters such as the distribution and interconnectedness of fluids and melts and hydrothermal alteration mineralogy (Muñoz, 2014; Trainor-Guitton et al., 2017). Existing MT studies at Aluto concentrated on regional aspects of the active rifting (Hübert et al., 2018; Dambly et al., 2022) or imaging and integrated petrological interpretation of the magmatic heat source and underlying crustal mush system (Samrock et al., 2020, 2021). This study presents the first detailed imaging of the geothermal system under Aluto (i.e., the depth range from the surface to the heat source).

The conceptual model for a magma-driven geothermal system is shown in Figure 2. The sketch builds upon earlier conceptual reservoir models, presented, for example, in Berktold (1983) and Pellerin et al. (1992). The model is explained as follows. Central to the formation of magma-driven geothermal systems, also known as high-temperature or high-enthalpy geothermal systems, is the existence of a powerful heat source (Jolie et al., 2021). In volcanic environments, the heat source consists of a shallow magma reservoir typically at depths of 3–6 km below the surface (km b.s.). During the course of fractional magma crystallization, volatiles exsolve from the melt phase. These hot and buoyant fluids rise along permeable structural damage zones, such as faults (e.g., Wallis et al., 2018). This can lead to rather complex temperature distributions with strong well-to-well gradients (Figure 2). The temperature field distribution will control the hydrothermal alteration mineralogy. Between temperatures of roughly 70°C–150°C, the dominant hydrothermal alteration product is smectite clay (Ussher et al., 2000). At higher temperatures, smectite becomes thermodynamically unstable and is gradually replaced by illite clays (Pytte and Reynolds, 1989). The hydrothermal alteration regime below 250°C is called argillic alteration, and the formed alteration products are typically referred to as clay caps. The bottom of the clay cap (⁠$T≤250°C$⁠) is commonly defined by the 0.1 S/m electrical conductivity isosurface (see, e.g., the conceptual reservoir model in Gíslason et al., 2015). At higher temperatures >250°C, the alteration regime becomes propylitic, and smectite clays are replaced by alteration products such as chlorites and epidotes. Smectite and illite clays are characterized by a high cation-exchange capacity, which makes the argillic alteration zone electrically highly conductive (Ussher et al., 2000; Lévy et al., 2018). Propylitic alteration products have low cation-exchange capacities and are therefore less electrically conductive. Hence, a temperature increase below the clay cap is correlated with a decrease in electrical conductivity.

Earlier studies, where the conceptual reservoir model (e.g., Pellerin et al., 1992) was adopted, focused on seeking out clay-cap-like structures and paid little attention to the electrical conductivity signature of the deeper magmatic heat source that drives hydrothermal convection and exerts controls on other parts of a geothermal system. However, this “clay cap hunting” alone was soon realized to be insufficient for identifying productive zones within a geothermal reservoir. Locating the magmatic heat source and mapping its electrical conductivity in 3D is essential to understanding the structure and the thermochemical state of a magma reservoir and to localizing upflow zones that form above magmatic intrusions (e.g., Bertrand et al., 2012; Samrock et al., 2018, 2021; Ishizu et al., 2022; Yamaya et al., 2022). Potential hydrothermal upflow zones can be identified by analyzing the clay cap topography (e.g., Anderson et al., 2000): because the temperature pattern will follow upflow zones, the temperature isosurfaces, and hence the argillic alteration front, will be deflected and pushed closer to the land surface, resulting in an updoming of the clay cap (see Figure 2). A similar idea was recently put forward by Trainor-Guitton et al. (2017) and we provide the first practical case study, which confirms its feasibility and can motivate its wider adoption in the future as a new characterization sign for volcano-hosted structurally controlled geothermal systems.

Phase tensors are often displayed in the form of ellipses (e.g., Booker, 2014), with the tensor’s singular values (⁠$Φmax,Φmin$⁠) defining the length of the major ellipse axes and hence the ellipticity. The ellipticity is a proxy for the strength of lateral conductivity gradients, and its direction indicates the geoelectric strike direction with a 90° ambiguity. The asymmetry of the phase tensor is quantified by the skew angle $β$⁠. For 1D and 2D subsurface electrical conductivity distribution, the phase tensor is symmetric and the skew angle is zero; for the general 3D case, the skew angle is nonzero (e.g., Caldwell et al., 2004). In Figure 3, we show phase tensor ellipses for two representative periods of the observed MT data on the Aluto. At short periods (Figure 3a), it can be seen that low skew values $β$ and symmetrical PT ellipses indicate low subsurface dimensionalities and increasing electrical conductivities with depth as suggested by high $Φmax$ values of above 45°. At longer periods (Figure 3b), higher ellipticities indicate that the geoelectric strike direction is clearly influenced by the intersecting AJFZ on the center of the volcanic complex (compare to Figure 1). High skew angles $β$ suggest that a 3D interpretation of the data is required.

For numerical 3D modeling and inversion of the MT data from the Aluto (Figure 1), we use GoFEM (Grayver, 2015; Grayver and Kolev, 2015; Arndt et al., 2020), which has been successfully used in previous 3D MT studies (e.g., Samrock et al., 2018; Käufl et al., 2020; Munch and Grayver, 2023). GoFEM is a 3D inversion code that is based on the finite-element technique and exploits adaptive mesh refinement techniques for forward and inverse modeling. Meshes created with GoFEM consist of nonconforming hexahedral cells, whose size accounts for MT site distribution and data resolution (Figure 4). Flexible meshing allows us to accurately incorporate real topography, which is crucial for modeling topography effects in MT data (Käufl et al., 2018; Soyer et al., 2019). We note that the terrain in the study area is characterized by a pronounced relief, with elevations that range from approximately 1600 to 2400 m above sea level (a.s.l.). Hence, an accurate representation of, and accounting for, topography during 3D modeling is crucial for resolving the structure of the geothermal system, most notably the clay cap, which is expected to occur in the upper 2–3 km b.s.

To obtain a 3D subsurface model, we inverted phase tensors. Prior to the inversion, we assigned a 5% error floor to the impedance tensor elements row-wise. The obtained uncertainty was subsequently propagated to the phase tensor. The mesh used in the inversion is shown in Figure 4. The starting model consists of a homogeneous half-space with a data-informed electrical resistivity of $ρ¯assq=26.3Ωm$⁠. Here, $ρ‾assq$ is the regional geometric mean resistivity over all sites and frequencies derived from equation 6. SSQ average of impedance tensors is known to be less influenced by galvanic distortions (Rung-Arunwan et al., 2017) and thus is used to guide the phase tensor inversion and ensure its robust convergence (Rung-Arunwan et al., 2022). This study differs from previous analyses of MT data from the Aluto volcano (Samrock et al., 2015; Cherkose and Mizunaga, 2018; Samrock et al., 2020) in that we included significantly more MT stations (178 in total, Figure 1) and added higher frequencies (the total range is 338–0.0058 Hz) to enhance the resolution in the shallower part of the model, where we expect the clay cap to occur.

In this section, we will first analyze the data fit of the model, then discuss the final conductivity model, and finally provide an interpretation of the structure signified by the conductivity. Finally, we will discuss the implications of the imaged conductivity structure for the geothermal system and prospective areas for future drilling.

The normalized root-mean-square (rms) misfit decreased from 2.88 to 0.93 within six nonlinear iterations, indicating that all data fit within the prescribed uncertainty. Histograms of the normalized residuals and crossplots of observed versus predicted data are shown in Figure 5. Note that the normalized phase tensor residuals for the initial model show a multimodal skewed distribution due to the nonlinear dependency on conductivity (Caldwell et al., 2004). As shown in Figure 5a, the residuals of the final inverse model exhibit a symmetric distribution with a zero mean. The crossplots of the observed and predicted phase tensor elements in Figure 5b show that dominant elements $Φ11,22$ follow the diagonal, indicating that responses from the final model have no systematic bias relative to the observed data.

Vertical slices through the obtained model are shown in Figure 6. A comparison with the conceptual reservoir model in Figure 2 shows that the key features of a magma-driven geothermal system are well recovered in the MT model. C3 is the magmatic heat source that was already imaged and described in detail in Samrock et al. (2021). The maximum melt content within C3 was calculated to be 10–15 vol%. The melt fraction was derived through an integrated MT and petrological analysis including thermodynamic modeling of the melt evolution and multiphase mixing laws considering a three-phase system consisting of solid, melt, and magmatic volatiles (see details in Samrock et al., 2021). In addition, a large amount of magmatic volatiles, exceeding 5 vol%, are suggested within C3. The very conductive layer C2 in the upper approximately 2 km b.s. under the central part of the volcanic complex is the prominent clay cap and the argillic hydrothermal alteration zone, for which well logs identified the presence of kaolinite, smectite, and illite clays, leading to a strong enhancement of the electrical conductivity (e.g., Teklemariam et al., 1996; Sisay, 2016). West and south in the flat plains surrounding the volcanic complex, high electrical conductivity values of >0.1 S/m are mapped in the thick layer occupying roughly the upper 2 km b.s. The fumarolic activity around well LA-2, west of the Aluto volcano (Figures 1 and 6a), and the clay-cap-like electrical conductivity signature led to the decision to drill well LA-2, which was found to be nonproductive. The high electrical conductivities under LA-2 are caused by a mix of electrically conductive lake sediments and montmorillonite clays (Hochstein et al., 2017). The maximum temperature in LA-2 is 105°C at 1.4 km depth and is related to advective thermal outflows that are underlain by colder rocks. The bottom of the well encountered impermeable colder basalt (note the related strong decrease in electrical conductivity at the bottom of well LA-2 in Figure 6a). The situation was similar for well LA-1, south of the Aluto volcano, which is also nonproductive. However, in both wells, high-temperature alteration minerals were found, which indicate much higher paleo-temperatures of up to 250°C. Hochstein et al. (2017) introduce the term paleooutflow zone to discriminate between outflows of different geologic ages that result from variations in convective upflows, caused most likely by alternating magmatic activity, volcanic eruptions, and magma recharge. To discriminate between the clay caps related to paleooutflows and active hydrothermal convection, it is therefore crucial to image the active magmatic heat source. The magma reservoir and heat source at the Aluto-Langano volcanic center occur at depths below approximately 3.5–4 km below sea level in the central part of the Aluto dome. Here, productive wells are all drilled into a narrow zone under the AJFZ, which was identified as the major zone of hydrothermal upflow (the red circle in Figure 7b).

Wells LA-3 to LA-10D, drilled atop the Aluto dome, show strong crosswell temperature gradients that indicate a rather complex flow pattern of hot circulating fluids. Although the AJFZ was identified as the major hydrothermal upflow zone, wells LA-7 and LA-5, drilled approximately 1 km west and east of the AJFZ, respectively, were not productive, with both wells encountering an outflow zone and even temperature reversals at depth (Woube, 1986; Benoit et al., 2007).

To understand the hydrothermal flow pattern and identify the upflow and lateral outflow zones, the clay cap topography needs to be mapped. The bottom of the clay cap experiences updoming over hydrothermal upflow zones and lowering over outflow zones (Trainor-Guitton et al., 2017), as the argillic alteration mineralogy is controlled by the temperature distribution. This effect can clearly be seen in Figure 7, which shows the depth to the bottom of the 0.07 and 0.1 S/m isosurfaces of conductor C2, identified as the clay cap. As shown, the clay cap experiences updoming and thinning along the AJFZ and clear clay cap lowering and thickening west of the AJFZ, where wells encountered the outflow zone and temperature reversals (see also Figure 8). Interestingly, another zone of clay cap updoming can be observed along the WAFZ west of the AJFZ (Figure 7). The reservoir below the WAFZ has not yet been considered for drilling, but the clay cap updoming revealed by the model suggests that it also might represent a productive zone because it is located along a potentially permeable geologic discontinuity above the magmatic intrusion. Fumarolic activity is also observed at the WAFZ (Figure 1), which provides further evidence that the WAFZ is a permeable upflow zone. For a detailed understanding of the inferred hydrothermal circulation pattern and conceptual reservoir model of the Aluto-Langano geothermal system, compare Figures 2 and 7b.

The presented study demonstrates the value added by high-resolution 3D MT imaging in geothermal exploration. Such imaging enables a geophysically guided characterization of magma-driven, high-temperature geothermal systems by providing key information about geothermal reservoirs, the depth and extent of the magmatic heat source, and the clay cap structure above the convective reservoir. Integrating this information with existing geologic and drilling data allows us to infer hydrothermal convection patterns, which eventually help in derisking the drilling. It is noteworthy that the conductivity model of the Aluto-Langano geothermal field shows why the first wells (namely LA-1 and LA-2) were not productive despite similar lithologies encountered in productive and nonproductive wells. The likely explanation is that hydrothermal alteration products encountered in nonproductive wells LA-1 and LA-2 formed in the past when higher paleotemperatures were locally present (Hochstein et al., 2017; Vidal et al., 2018). Therefore, to understand the present-day circulation of geothermal fluids and to optimize the selection of geothermal well locations, it is crucial to consider hydrothermal alteration products in connection with a magmatic heat source. Indeed, no evidence for the presence of notable amounts of melt below wells LA-1 and LA-2 can be found in our model (Figure 6). The model reveals only low electrical conductivities below these nonproductive wells (Figure 6). These low conductivities are indicative of cooled crystalline rock, which was encountered in logs at the bottom of the wells, with bottom-hole temperatures below 105°C (Hochstein et al., 2017). We note that LA-1 and LA-2 are at the edge of the main survey area, yet the data still provide sensitivity to the conductivity structure underneath these wells. We performed numerous modeling experiments (not shown here) whereby the final inversion model was modified underneath these regions and the rms values per each site were calculated for the modified model and compared against the rms values for the inversion model. These tests confirmed the significant sensitivity of the MT data to artificially added conductors (mimicking a heat source) below LA-1 and LA-2. Furthermore, a previous study by Samrock et al. (2015) had an MT site next to LA-1 but did not find signatures of melt.

Productive upflow zones typically form above the intrusion within structural damage zones, such as faults, and can be identified by clay cap thinning and upwarping. Note, in studies about the value of geophysical data for geothermal exploration (Trainor-Guitton et al., 2017; Trainor-Guitton, 2020) this correlation was identified at other prospects as a potentially important characteristic for locating high steamflow (i.e., upflow) regions in geothermal reservoirs and is confirmed here. This study thus paves the way for a future avenue of research that can further strengthen and refine this new exploration metric for prospects with geologically similar settings. We further observe that conductor C1, below the LA-1 and LA-2 wells, is significantly thicker and exhibits a more homogeneous electrical conductivity distribution with a flat bottom, unlike conductor C2 under the central part of the Aluto volcano, which clearly shows updoming and downwelling. This observation is likely the consequence of downward propagation of the argillic mineral alteration front that occurred during the cooling of the paleogeothermal reservoir.

MT exploration for magma-driven geothermal systems offers a high value of information by derisking exploration well drilling provided the following conditions are met. First, mapping the heat source below the clay cap is key to distinguishing between nonproductive and productive zones in a geothermal reservoir. Second, we observe a clear correlation between upflow zones and geologic discontinuities. The upflow zones can be inferred by mapping the regions of clay cap thinning and updoming. To achieve the aforementioned objectives and facilitate the successful siting of productive wells, physically consistent 3D conductivity modeling of MT data and a dense MT site distribution that covers the entire extent of a geothermal prospect are of utmost importance. The case study of the Aluto volcano is likely relevant not only for other volcano-hosted geothermal sites within the MER but also should be tested for other high-temperature prospects, many of which expose similarities in geologic settings.

We thank the editor E. Gasperikova, reviewers G. Hill, K. Christopherson, M. Comeau, and an anonymous reviewer for their comments that helped us improve the original manuscript. We thank S. Fisseha, G. Endarge, and B. Cherkose for their collaboration and strong support. Furthermore, we thank E. Layman and W. Cumming for valuable discussions about the 3D model. F. Samrock thanks A. Jackson and A. Kuvshinov for initiating the ETH MT project at the Aluto. M. O. Saar and F. Samrock thank the Werner Siemens Foundation for their endowment of the Geothermal Energy and Geofluids group at ETH Zurich. A. Grayver was supported by the Heisenberg grant from the German Research Foundation, Deutsche Forschungsgemeinschaft (project no. 465486300). M. L. T. Dambly was supported by ETH grant ETH-02 19-1. The 3D inversions and modeling were carried out at the Swiss National Supercomputing Center under project ID s1106.

The MT data collected at Aluto by ETH Zurich are available via the IRISEMTF Database (http://ds.iris.edu/spud/emtf) under the Project ID Ethiopia.R1.2012. The MT data by the Geological Survey of Ethiopia are available for academic purposes on request from the Geological Survey of Ethiopia, as was the case for this study.

Biographies and photographs of the authors are not available.