Geothermal Mechanism and Reservoir Characterization with Multiscale Analysis of the Potential Field and MT Profile in the Yanqing-Huailai Basin, North China

The city of Zhangjiakou, located in north China and northwest of Beijing, will be one of the host cities of the Olympic Winter Games 2022. The exploration and exploitation of geothermal energy are one of the necessary approaches to low-carbon and green Olympics. The study area covers the middle part of Zhangjiakou city geographically, and in the division of the tectonic unit, the Yanqing-Huailai basin (YHB) is included. In the basin, abundant geothermal resources were discovered and field survey was carried out for mechanism research of geothermal anomaly and further exploitation.

According to the survey and research results, tens of geothermal anomalies are discovered and their distribution is concentrated around villages like Houhaoyao, Wanyao, and Xijiabao. The hot springs and geothermal boreholes have a medium and low temperature from 25°C to 89°C, and the corresponding depth is usually less than 500 meters [1–3]. In addition, the YHB situates in the structural convergence position of the northeastern part of the Fenhe-Weihe graben and the northwestern part of the Zhangjiakou-Penglai fault zone (ZPF) [4], which causes fault activities and frequent earthquakes. The development and occurrence of earthquakes might be related to the distribution of geothermal manifestation [5–7]. In the process of infiltration and deep circulation of atmospheric precipitation and surface water, the temperature will increase due to the influence of geothermal warming and residual heat of magma. The rock is broken, and the tensional pinnate structure is developed at the intersection of faults, which forms the enrichment zone and upwelling channel of geothermal fluid [3, 8].

The field survey work aiming at geothermal energy and earthquake occurrence in the YHB during passed decades has accumulated abundant data and materials. However, the lack of intensive study and understanding of geothermal anomaly restricted the development and utilization in the region. It is recognized that the reservoir in YHB includes Archaean gneiss, Jixianian dolomite, and Yanshanian granite principally but the spatial distribution of reservoirs and information in a deeper position have formed yet. The formation mechanism and key controlling factors remain ambiguous [1–3, 9]. The undulation of the Curie isothermal surface, as a result of the combined effect of convection, local radioactive heating effect, and upper mantle thermal radiation, is closely related to deep thermal events [10, 11]. It is necessary to further shed light on the indicating the significance of the Curie surface depth (CSD), which generally was calculated from regional magnetic data with susceptibility interface inversion approaches [12, 13]. Besides, magnetotelluric (MT) and seismic profiles were frequently used for geothermal and earthquake study in the YHB area to obtain deep information [2, 14, 15] but the lack of potential field application and integrated interpretation could not satisfy the need of comprehensive understanding and three-dimensional modelling of a regional geothermal reservoir [16–18].

In this paper, regional gravity and aeromagnetic data are interpreted with wavelet multiscale analysis and fast imaging methods. The resistivity inversion cross-section of MT data, CSD data, geothermal gradient, and terrestrial heat flow statistics is collected and involved in the integrated interpretation for the geothermal mechanism and three-dimensional model. With the advantage of high vertical resolution, MT profile data is primarily used to depict the characteristics of strata and properties of faults. Anomaly separation of potential field data is an essential step before the inversion and interpretation of geological targets. Data observed on a measurement plane with various altitudes have different definitions and in result reflect the characteristics of geological targets with coarse to fine scales and deep to shallow depth [19, 20]. In our study, the potential field data are decomposed into multilayers, and on each layer, the horizontal distribution of faults, reservoir, and intrusive bodies is delineated. The depth of the source that causes the anomaly would be estimated with a regularized downward continuation method [21], and in the process, resistivity inversion cross-section provides a priori information constraints. The data of temperature and heat flow play a key role in analyzing regional geothermal geological conditions and determining the geothermal objective area [11, 22, 23].

The scheme of integrated interpretation with multisource and multiscale geophysical data is proposed in our study. On the basis, the deep genesis of the geothermal anomaly would be revealed and the local geothermal reservoir in the YHB is characterized. At last, a conceptual geothermal mechanism model is built for a better understanding.

The study area is located in the overlapping compound position of the western part of the Yanshan intraplate orogenic belt and the northeastern boundary of the Fenhe-Weihe Graben (FWG) system, the northern part of the Sino-Korean quasi platform [4]. The Yanshan intraplate orogenic belt experienced intense structural deformation, magmatic and volcanic activity, and synorogenic sedimentation during the Jurassic-Cretaceous Yanshan Movement (Figure 1(a)). In addition, the FWG system has gradually formed since the Tertiary and it appears as a strongly active seismic zone. The upper mantle material in the graben system is uplifted, and the crust tends to be stretched. The development of the graben system controls the formation of the YHB, and at the same time, the structural pattern of the study area before the Cenozoic has been greatly modified [24]. The tensional environment and active Cenozoic tectonic movement are attributed to the formation of regional geothermal manifestations.

According to the regional geological survey report and boreholes in the surrounding mountainous areas and the basin, all strata from the Archean to the Cenozoic were discovered except for the Lower Paleozoic, Upper Ordovician, Upper Paleozoic Silurian, Devonian, and Carboniferous-Permian [25]. Most part of the study area is covered by Quaternary strata. The sedimentary thickness is generally tens of meters and the thick part could be up to 150–700 meters [3]. Neogene strata was developed and complete, including Miocene and Pliocene. Sandy conglomerate and sandpaper clay are interbedded, with no exposure on the surface. The Cretaceous system is absent in this area, and the Jurassic is widely distributed throughout the region, with a main lithology of sandstone, conglomerate, tuff, and andesite. The exposure of the Triassic strata is not widespread. Archaean is mainly distributed in the Pangjiapu, Mayukou, Tumu, Houhaoyao areas, etc. It is a set of metamorphic rock series formed by mid-deep-deep regional metamorphism [26].

Following the formation of the Mesozoic NNE-NE trending main fault system in the vast areas of eastern China, a large number of NW-NWW compressional (or compression-shear) and associated tensional and tension-shear fault systems have developed since the Cenozoic. It forms the NW-NWW faulted tectonic system [27–29]. Among them, the Zhangjiakou-Penglai fault zone (ZPF) penetrates the Cenozoic Bohai Bay Basin, the North China Craton basement, and the Mesozoic Yanshan Orogenic Belt along the NWW direction and has played a key role in controlling the deposition, magmatic activity, and tectonic development of different evolution stages on both sides of it [30]. Frequent earthquakes in the recent years indicate that ZPF remains seismically active, and focal mechanism solutions suggest that its slip vector is sinistral [31–33]. The main faults and the information of their property and activity are listed in Table 1. Under this structural background, the Late Cenozoic in the Zhangjiakou area underwent multistage extensional structural deformation and developed a series of intermountain fault basins, such as the Yanqing-Fanshan Basin (YFB), Huailai-Zhuolu Basin (HZB), and Zhangjiakou-Xuanhua Basin (ZXB).

Geological structures have an obvious controlling effect on the distribution of geothermal manifestations (hot springs and geothermal boreholes) according to our survey statistics. For instance, the Houhaoyao geothermal field is bounded by five faults [3]. The NW direction fault lying in the middle part plays a major role in thermal control, and to its southwest side, an intersection of faults along the NEE and NNE directions causes rock disintegration and development of plume-like fractures. The fracture zone is advantageous to the formation of geothermal water gathering and upwelling.

The study area is part of the Daxinganling-Taihangshan structural magma belt and the magmatic activity is strong. The main activity era is Mesozoic, and the strong period happened in the middle and late Mesozoic. The magmatic rock mass formed during this period has the characteristics of large scale, wide distribution, diverse shapes, and significant mineralization. It was mainly produced in the shape of rock foundation, rock strain, rock wall, rock branch, lava blanket, and rock sheet with different scales. The magma intrusion and volcanic eruption are controlled by the Dahenan-Chicheng deep fault (F₄ in Figure 1(b)) spreading from NE to NNE.

The study area belongs to the intermountain basin and piedmont zone, where an insular or half-insular fracture zone geothermal reservoir with linear and banding distribution is developed. The types of groundwater in the study area are divided into interstitial water in loose rocks, karst water, and bedrock fracture water according to water-bearing media. Interstitial water distributes along the lower reaches of Sanggan River, where advantageous water-proof layers are scarcely produced and gravel is exposed on the surface [3]. The aquifers are mostly an exoteric single-layer structure with strong permeability. Precipitate water and surface runoff vertically infiltrate into the ground through sand and gravel, and mix the water from mountains along the basin margin and river to recharge the runoff area for groundwater. Geothermal water circulates up generally through deep faults and occurs in fault zones and fracture developing zones. The upper part of the bedrock reservoir is covered by the Quaternary stratum, which constitutes the caprock of the bedrock reservoir and has the effect of heat accumulation and preservation for geothermal water. Meanwhile, the groundwater in the Quaternary sand and gravel layer is heated by thermal conduction to form the Quaternary geothermal reservoir.

The geothermal gradient distribution in the study area is controlled by the geological structure pattern of the unevenness. The relatively high temperature zone is consistent with the structural uplift zone. The caprock geothermal gradient in the Houhaoyao area is greater than 3.5°C/100 m and could reach 4.0~5.0°C/100 m or more locally [3]. Vertically, it is mainly controlled by the lithology of the formation and the geothermal gradient of rocks with high thermal conductivity is smaller. Therefore, the geothermal gradient of the Cenozoic is generally greater than that of the underlying bedrock. Terrestrial heat flow is a comprehensive thermal parameter, which can reflect the thermal state of a region than basic geothermal parameters (such as temperature and geothermal gradient) more accurately. With regard to the study area, it has the characteristics of a higher value in the transition zone between the mountain and the plain and the Cenozoic graben than in the mountain area [7, 22]. The statistical data are listed in Table 2.

Geophysical properties of rocks in the study area are measured and collected (Table 3, figures in the bracket denote the number of rock samples). The statistics play a role of geophysical interpretation foundation, and from Table 3, it could be discovered that distinct differences in physical properties between different strata exist. It is worth noting that the crust of the YFB and HZB is a commonly high- and low-velocity alternate structure according to seismic imaging. The middle-upper crust undergoes velocity reversal frequently, and there is a high-dip fault zone extending to the Moho surface. This may be regarded as a deep sign of the gestation of earthquake geothermal anomalies in this area.

Satellite gravity data is from the new global marine gravity model from CryoSat-2 and Jason-1, that is, S&S Global Anomaly V29.1 downloaded at https://topex.ucsd.edu/cgi-bin/get_data.cgi. The gridding interval is $1′×1′$⁠. Regional gravity data has a resolution of $2km×1~2km$ and provides the information of strata and fault distribution in the three-dimensional space, while aeromagnetic anomaly is mainly used to demonstrate the characteristics of intrusive rock mass, which was speculated to contribute to the local geothermal anomaly.

The aeromagnetic data covers an area of Xuanhua and Huailai with grid point intervals of $0.5km×0.5km$⁠. The data were collected by a fixed-wing aircraft flown at an average altitude of 500 meters, and the total anomaly has been reduced to a pole. The magnetic measurement program was carried out by the Hebei Institute of Geological Survey in 2016 aiming at exploring and assessing the mineral resources. Linear geological structures and high- and low-anomaly zones are demonstrated with the total field data to interpret prospective provinces and mineral plays.

The MT profile was acquired and interpreted by Jilin University and has been published [2]. On their basis, we modify and reinterpret the profile to obtain a more explicit understanding on the resistivity structure and the controlling effect of low-resistivity anomaly over the formation of regional geothermal energy. In our field survey, the information of the temperature, depth to water table, strata of hot springs, and geothermal wells was collected. In addition, physical properties of rocks, geothermal gradient, and terrestrial heat flow data [3, 7, 22] were gathered into the integrated interpretation.

Since the temperature of groundwater in the study area varies from 12.4°C to 82.0°C, its spatial distribution shows a decreasing trend centered on the Houhaoyao and Xijiabao geothermal field. In order to better distinguish the geochemical characteristics of geothermal fluids, all samples are classified by temperature as high temperature (>45°C), medium temperature (25–45°C), and low temperature (<25°C). We analyze the hydrogeochemical parameters of the samples according to the abovementioned classification. The range of high-temperature geothermal water samples is 45–82°C, and the average value is 54.53°C. Sampling sites are dominantly located in the central area of the Houhaoyao and Xijiabao field. Testing parameters are counted and summarized in Table 4.

The nonuniform density distribution exists in the crust and upper mantle, which is related not only to the composition but also to the internal temperature variations in the lithosphere. The density of the lithosphere in north China is significantly nonuniform in both the horizontal and vertical directions and exhibits a notable segment-wise spatial distribution pattern. Two tectonic units, namely, the Taihangshan tectonic zone and the FWG, constitute a gravity gradient transition zone striking the NNE direction, which traverses the YHB and produces a significant difference in the density distribution on the two sides [35]. On this basis, satellite gravity anomaly provides an approach and window to analyze the density inhomogeneity in the lithosphere. It is obvious to discover high values along the NE direction and low values along the NW direction in the gravity map (Figure 2); the trends were considered as a combined effect of FWG and ZPF. Besides, the study area lies at the margin of Yan Mountain and North China plane to the southeast. The variations of anomaly correspond to the gravity gradient transition zone, where deep and major faults develop and the crustal movements are usually intensive. The conditions would provide a favorable geothermal background.

In this paper, the structure of the sedimentary basin refers to the transverse density and magnetic susceptibility inhomogeneity of different depth levels and the scales of potential field anomalies are dependent upon sizes of the sources located in different depths which vary accordingly [39]. Based on the geological characteristics of the YHB, we inversed the potential field perturbation corresponding to the three sets of strata that covered the major part of the basin. Daubechies wavelet is applied into the data decomposition and reconstruction with the advantages of orthogonality, compact support, and approximate symmetry [38]. Gravity and magnetic anomalies on each scale represent the response of sources at different depths as inverse relation exists between the scale and the frequency. It is not necessary to set the range of depth related to anomalies at each scale, and besides, the interpretation and inversion are applied without depth weighting. The method could improve the resolution in the depth direction effectively and provide more detailed deep information.

Furthermore, two regional gravity profiles are drawn in Figure 2 to study the vertical distribution of strata and faults using 2.5D inversion. It should be noted that the profile data is extracted from regional gravity instead of satellite gravity. Figure 3(a) presents a profile striking near the SN direction and passes through a couple of faults striking the NE, NNE, and NW directions. It is obvious that the Quaternary Fanshan Basin and Zhuolu Basin are controlled by normal deep faults at the south and north ends. At the south end of the Fanshan Basin, a rock mass developed and it is believed to form during the period of the Yanshan Movement. Figure 3(b) shows the inversion result of the profile BB$′$⁠. The profile strikes 45°, and the Taohua-Nuanquan Basin is a Cenozoic monoclinic basin controlled by the Yuxian-Yanqing Fault, whose northwestern part descended and deposited a huge thickness of Quaternary sediment. The Zhuolu Quaternary Rift Basin developed in the central part and is a Cenozoic faulted basin, and the northwest side is controlled by secondary normal faults. The southwestern end is twisted by the Xiahuayuan Major Fault and is separated from the southern depression.

Stable downward continuation is a promising technique for estimating the depth of sources. Based on the filter curve in the wave number domain, we developed a new depth estimation method using the Chebyshev-Padé approximation function. The method was published [21, 40] and utilized successfully in the interpretation of potential field data.

The potential field data are decomposed and reconstructed on four layers with the wavelet multiscale analysis method introduced in Section 3.2. The anomalous field on each layer denotes the distribution of density or magnetic inhomogeneity at different depths. Figures 4(a) and 4(b) represent the results of aeromagnetic and Bouguer gravity anomaly, respectively, and from them, the information of the thickness of strata, penetration depth of faults, and the scale of rock mass could be extracted. What is mentionable is that gravity anomaly reveals more detailed information than magnetic data on most layers in the study area. Gravity data contains evidently more features associated with tectonic units, faults, and strata (Figure 3), while magnetic data reflects the spatial distribution of magmatic rock and metamorphic rock more clearly according to Table 2. Additionally, the temperature measurements of field survey points are superimposed, which provides an approach to evaluate the relations between geothermal anomaly and potential field features. In order to demonstrate the residual density inhomogeneity of the upper crust in the three-dimensional space, a fast-imaging method based on downward continuation of potential field data is applied [21]. A set of three-dimensional data visualization patterns, including face render, contour line, and contour surface, is utilized, and the result is shown in Figure 5.

Magnetotelluric sounding data was collected along the profile starts from nearby the Houhaoyao geothermal field to the northwest and has a length of 80 km. Figure 6(a) is the apparent resistivity inversion result, and the interpretation depth reaches the asthenosphere. Based on the geophysical statistics in Table 2, the apparently layered resistivity structure could be identified. The first layer, above 2 km, corresponds to the sedimentary cover, which functions as thermal insulation. The high-resistivity layer is inferred to the basement composed of carbonate deposits and Archean gneiss. The third layer, under 25 km, is mainly characterized by a high conductor extending to a depth of 100 km, with a resistivity value of 3~30 Ω·m [2]. The inversion result of the MT profile shows that the high-conductivity anomaly has a fair deep root, which penetrates the lithospheric mantle and connects with the asthenosphere. There may be an upwelling of mantle-derived material in this area, which could be one of the controlled factors of regional geothermal anomaly [15].

In Figure 6(b), it is clear to see that the uplift of CSD spreads a majority part of the YHB and the peak appears near Zhuolu. Three main factors that affect the depth in the inner crust are roughly as follows: local geothermal gradient, the composition of the titanomagnetite in the crust, and the internal temperature changes with the hydrostatic pressure [11]. Among the abovementioned factors, the temperature in the crust has the greatest influence. Therefore, the undulation of CSD reflects the distribution of the temperature field in the deep crust and it appears to be the shallowest interface in analyzing the thermal structure of the lithosphere. The CSD in the study area also shows a block-like structure between the uplift and the depression, with a depth range of 19~29 km. Although its undulation is affected by thermal conduction of the Moho interface, it is more importantly controlled by the thermal convection of the crustal fracture. In terms of morphology, it does not correspond to the geological tectonic unit each to each.

Figure 6(c) shows the CSD along the MT profile, and it is obvious to see that the trend is in accord with thermal upwelling. It is demonstrated that the survey points with high temperature around the Houhaoyao and Xijiabao area correspond to a significant gradient zone of CSD, where obvious geothermal anomaly, stress reduction, and energy release are discovered. In MT cross-section, it is also clear to see the low-resistivity thermal upwelling migrates to the shallower part of the crust. Under the combined action of this thermal upwelling and the deep faults, geothermal anomaly could be tracked.

In Figure 7, the piper diagram shows that the hydrochemical type presents various distributions. The hydrochemical characteristics of surface water show that the ion content is similar, Mg²⁺ is much higher than that of groundwater, and the hydrochemical type is Na•Mg-Cl•HCO₃•SO₄. Compared with reservoir samples, river water samples are closer to groundwater in the piper diagram, indicating that its chemical characteristics may be affected by groundwater recharge to some extent. The hydrochemical analysis of the temperature-based groups, the distribution characteristics of TDS, and various ion contents can infer that the geothermal water undergoes deep circulation and is mixed with shallower cold water.

From the diagram in Figure 8, cold groundwater samples are close to dolomite and high-temperature water samples approach anhydrite. Additionally, the ratio of Ca/SO₄ in the high-temperature samples is less than 1 and it is also smaller than those of the low-temperature and medium-temperature groundwater. Combined with the significantly low content of Mg²⁺ in high-temperature geothermal water, it is speculated that Ca ion concentration is greatly reduced due to the precipitation of minerals such as carbonate and dolomite. The conclusion is consistent with the observed carbonate precipitation and alteration of geothermal water that is prevalent.

The results in Figure 4 indicate that the deep and shallow subanomalies reveal the relation with ancient fluid active zones. The negative density disturbances inversed from the shallow subanomaly are mainly caused by Mesozoic fluid active zones, whereas the negative density disturbances form the deep subanomalies are mainly correlated with Paleozoic fluid active zones. The position of Dahaituo rock mass corresponds to low gravity and high magnetic anomaly. In Figure 4(a), the geothermal anomaly is distributed at the ends of the intrusive rock mass and dyke, so the formation is inferred to be related to magmatic activity. In Figure 4(b), the geothermal manifestation with higher temperature shows that it is commonly distributed along the gravity gradient belt, showing the correlation with basin margin faults. Combining gravity and magnetic information, it would be found that the distribution of geothermal anomaly shows obvious regularity.

At the four layers in Figure 4(b), the geothermal anomaly corresponds to gravity values from −100~−110 mGal, which denote the underground source with a specific residual density and could be interpreted as geothermal reservoir qualitatively. To further demonstrate the three-dimensional characteristics of this geological target, the 3D fast imaging of gravity data using regularized downward continuation is calculated and shown in Figure 5. The regularity in Figure 4 and the spatial feature of geothermal reservoir are integrated to build a 3D geothermal geological model.

In Figure 6(a), the resistivity feature of interest is the high-conductivity body and it was indicated to be an upwelling of the mantle-derived material in this area [41, 42]. The electrical structure shows that the heat carried by the thermal upwelling would be transferred to the surface passing through the high-resistivity basement. Accordingly, the deep faults and their secondary faults and rock fracture zone could provide the migration channel. So, the two aspects constitute the main controlling factors of geothermal manifestation.

For the CSD in Figure 6(b), Houhaoyao and Xijiabao geothermal fields are both located along the CSD gradient belt (Figure 9(b)). At a depth of more than 10 km in the crust, blocks on both sides of the belt may have a temperature difference of 100 to 200°C [11]. The produced thermal stress might drive deep fluid to migrate from the deep part to the shallower along the crustal major fault. The mechanism is in accord with the electrical structure.

Based on the abovementioned analysis, the temperature survey points and temperature distribution contours in the study area are combined with the CSD and geological information such as strata, fractures, and rock masses obtained from geophysical data is integrated to estimate the depth of the geothermal reservoir (Figure 9(c)).

According to the hydrochemical characteristics of the geothermal field, the water is dominated by a SO₄-Na type and is rich in fluorine (F⁻) and soluble silicate (H₂SiO₃). The presence of SO₄²⁻ ions is formed by the dissolved water of H₂S in the deep part, and the presence of F⁻ is caused by the influence of magma or hydrothermal fluid on hot water. Soluble metasilicic acid is generally considered to be related to the migration of SiO₂ in magma. The abovementioned information indicates that the formation of underground geothermal water is related to the hydrothermal activity caused by the deep magma source.

According to the differences in structural causes and thermal transfer pattern, hydrothermal geothermal resources can be divided into two types: sedimentary basin conduction and uplift-mountain convection. Conductive geothermal resources in sedimentary basins are mostly developed in plains and intermountain basins. The heat in the deep crust is continuously transferred upwards mainly through the conduction of the stratum. Convection-type geothermal resources in the uplifted mountain area take convection as the main heat conduction pattern, and a large number of developed faults in the mountain area create advantageous conditions for hot water convection. Our study shows that the geothermal type in this area could be uplift-mountain convection typically.

The Houhaoyao geothermal field is a typical one developed in the graben basin. It is located on the northwest side of the NE Dahenan-Chicheng Fault (F₄) and the NNE Yuxian-Yanqing Fault (F₅) and the north side of an anticline. The geological structure in the geothermal field is complex, and most of the faults are concealed. The NW trending fault in the central part is the main thermal-controlling fault. Due to the influence of NE trending and NNE trending faults, rocks are broken at the intersection and tensional plume-like structural fissures are developed, forming a geothermal water enrichment zone and upwelling channel. The distribution of the high-temperature zone is consistent with fault strike [3]. The strata from old to new are Archean gneiss, Proterozoic Changchengian system marine sediments, Mesozoic Jurassic sintered tuff, and Cenozoic thick sediments.

The upper reservoir of geothermal fields in the study area is shallow and acts as a channel for the release of deep thermal energy, which is formed by the conduction of hot water from the deep part to the medium. The previous speculations about the geothermal source of underground hot water include residual heat of magmatic rock, heat of geothermal gradient accumulation, heat of radioactive element disintegration, and mechanical heat of neotectonic movement [43]. The test results of the radioactive heat generation rate (RHGR) of rocks show that RHGR of gneiss, granite, and porphyritic granite in the study area is 0.694 μW/m³, 5.459 μW/m³, and 5.785 μW/m³, respectively. The heat source in the study area should not be dominated by the residual heat of the magmatic rocks, but by the heat conduction and temperature increase of the earth.

According to the geothermal geological simulation results of the Houhaoyao-Xijiabao area, referring to the regional geothermal geological background and existing self-heating geological data, the geothermal origin mechanism map of the Houhaoyao-Xijiabao area is drawn, as shown in Figure 10. As can be seen, the Houhaoyao geothermal field is located in the Yuxian-Yanqing Deep Fault (F₁) zone. According to the results of geophysical exploration in 2019, there is a north-northwest-trending fault (F₂) passing through the geothermal field. From the perspective of its formation and hot water movement mechanism, a shallow thermal storage is a shallow pore-type thermal storage formed by hot water flowing from the deep thermal storage to the Cenozoic stratum along deep faults and mixing with cold water for heat exchange. A deep heat storage is the basis and prerequisite of the shallow heat storage and is the source of the shallow heat storage. A shallow heat storage is the heat energy release channel and hot water discharge channel of the deep heat storage. The geothermal fluid is rich in metasilicic acid and fluorine and has a high salinity, indicating that the geothermal water dissolves a large amount of rock and soil components in the deep high-temperature and strong-pressure environment. Compared to those of Xijiabao, the metasilicic acid and fluorine contents of the Houhaoyao geothermal field are higher. Accordingly, the fluid circulation depth is deeper and the possible rock mass lies in the environment with higher temperature and stronger pressure.

Through the analysis of the geothermal geological model in the key research area, it is supposed that the Houhaoyao-Xijiabao geothermal field is formed under a ternary thermal accumulation mechanism. One is the conduction of deep heat flow, and the second one is the accumulation of thermal refraction effect on the shallow uplift area. The last one is the accumulation of groundwater convective heat, including convective accumulation of the thermal conduction fracture or contact belt of rock mass and groundwater migration. To be specific, the atmospheric precipitation in the Houhaoyao-Xijiabao area infiltrates into the ground from Lang Mountain and Laojun Mountain. In the process of deep circulation, it accepts the conduction and thermal accumulation of the deep heat flow and then flows through the thermal conduction fracture zone or the rock contact zone. In the bedrock fissure reservoir, it receives convection thermal accumulation during the migration process. Part of the hot water upwells and connects with the shallow aquifer and mixes with cold water for heat exchange, forming a shallow pore-type reservoir. In the raised area of the bedrock, it also accepts the heat refraction effect due to the high thermal conductivity.

The overall objective of our survey and study is to make full use of regional geothermal energy. The total amount of exploitable geothermal resources in the study area is abundant, and the flow rate is considerable. Water temperature ranges from 25°C to 90°C. However, the conventional exploitation way of irrigation and bathing might cause waste and pollution. On the basis of the geothermal accumulation mechanism and distribution characteristics, geothermal resources with higher quality would be discovered hopefully. The knowledge of the geothermal system could attribute to the development of a cascade exploitation route with high efficiency and sustainability.

Hydrothermal settings are common in tectonically active continental regions, where active faulting and fracturing provide possible pathways for the circulation and convection of geothermal fluids. The deposition and alteration of geothermal fluids and the newly formed mineral associations are generally considered to be produced under the action of medium- and low-temperature hydrothermal fluids, indicating that the thermal fields in the region are related to the residual heat of deep rock mass or thermal upwelling. The geothermal reservoir model recovered from potential field data and MT profile reveals the three-dimensional image of a hydrothermal system in the YHB. The model delineates reservoir characteristics and fault-controlled circulation of geothermal fluids. It offers new insights into migration processes of geothermal fluids in graben and fault zones. The results of this paper furthermore represent an improved mapping of the convecting geothermal system that forms under the impact of both neotectonic movement and magmatic intrusions. Exploitation and utilization of this renewable energy resource is currently increasing and could contribute to a low-carbon sustainable energy mix for the fast-growing green economies and carbon reduction strategy of China.

The data are not publicly available due to regulations of the related programs.

The authors declare that they have no conflicts of interest.

We wish to thank Dr. Liu Wenyu for providing the magnetotelluric profile and Curie isothermal surface data and discussing the scheme of integrated geophysical interpretation with regard to the geothermal mechanism. We also thank Dr. Qi Bangshen for the borehole materials and analysis on the relation between geothermal anomaly and earthquake development. This work was supported by the Geological Survey Project (grant numbers DD20190129 and DD20221677) and the Basic Research Funds (grant number JKY202007).