There are 1,000 geothermal wells with depths of 1500–2000 m drilled into the upper Tertiary sandstone geothermal reservoir in the Guantao Formation in the North Shandong Plain, which support about 61 million square meters of domestic space heating. Concerning the large-scale exploitation of the Guantao Formation geothermal reservoir and its important role in clean space heating, it is essential to fully understand the characteristics of the geotemperature, geothermal water flow, and hydrochemical features, especially the heat accumulation and geothermal water enrichment mechanism. Based on the data from pumping tests, temperature logging, water level measurements, and hydrochemistry and isotope sampling and analysis, in this study, the heat sources and accumulation, water sources, and enrichment mechanism of the sandstone geothermal reservoir in the Guantao Formation were explored. The geothermal gradient of the cap rocks and the reservoir’s temperature revealed that the geothermal temperature fields are controlled by fault structures, and the relative high temperature zones are distributed in the area where the uplifts and deep heat-controlled faults converge. The geothermal water-rich area is located in the ancient river channel zone. The closer the water-rich area is to the ancient river channel belt, the deeper and coarser the buried conglomerate aquifer is, and the better the water-rich condition is. The best water-rich areas are found in the alluvial fan of the ancient channel belt.

As a renewable and clean energy source that integrates heat and water resources, geothermal resources have received increasing attention due to the current energy shortage and the increasing awareness of environmental protection [16]. Although geothermal energy is a renewable energy resource, it is also a limited resource, and its recharge process is extremely slow. Many studies have focused on the geothermal potential. Therefore, it is of great significance to determine the thermal storage conditions, the characteristics of the water and heat sources, and their enrichment mechanism for scientific planning and sustainable exploitation of geothermal resources. Numerous scholars have also carried out research on these topics to achieve these goals [712]. For instance, a multidisciplinary approach was used to assess the geothermal potential of the Tønder area, North German Basin [13]. The permeability of the low-enthalpy geothermal reservoirs in the upper Triassic–Lower Jurassic Gassum Formation in the Norwegian-Danish Basin was predicted [14]. It is even more important to consider how the permeability of sandstones is affected during diagenesis since the fluid flow in reservoirs mainly depends on the permeability [1518]. The chemistry and isotopic compositions of water represent a large energy potential in the low-enthalpy geothermal system in the Cimino-Vico Volcanic District, Italy [19].

The purpose of this study was to identify the heat energy and geothermal water accumulation areas based on the characteristics of the geotemperature field, the hydrochemical field, and the geothermal water flow field. These accumulation areas are the target for geothermal energy development.

As a typical reservoir in the North China Plain, the Guantao Formation sandstone/glutenite reservoir in the North Shandong Plain represents a valuable case study (Figure 1). There are several different tectonic units in the North Shandong Plain, which have different temperatures, total dissolved solids (TDS), special chemical compositions, and water abundances. Among them, the Guantao Formation geothermal reservoir contains the most abundant resources and has the highest degree of exploitation and utilization. Correctly, understanding the characteristics of the geothermal field and the hydrochemical field and the enrichment mechanism of the geothermal water is of great significance to the scientific exploration, well determination, rational reinjection, sustainable development, and utilization of the geothermal resources in the future [8, 1012, 2023].

2.1. Tectonic Setting

The study area is the North Shandong Plain, which is located in the southern part of the Bohai Bay Basin (Figure 1). The Bohai Bay Basin is located in the northeastern part of northern China. It is bounded by the Tancheng-Lujiang fracture zone to the east, the Taihang Mountain eastern fracture zone to the west, the Yanshan orogen to the north, and the western Shandong uplift to the south, with a total area of about 200,000 km2 [24].

The North Shandong Plain is a Mesozoic-Cenozoic fracture-depressed basin developed in the North China Platform. In terms of tectonic units, the plain is located in the North China Depression in northern China Plate. Since the Mesozoic-Cenozoic, under the influences of the Himalayan Orogeny and Yanshan movements, the fracture structure was developed, and the concave-convex tectonic units were formed, namely, the Jiyang depression and the Linqing depression (Figure 2) [2527]. Affected by the movement of northern China in the Neogene, the rifting activity basically stopped, and the basin gradually entered the subsidence stage. The Late Quaternary Himalayan movement caused the North China Depression to rise as a whole and to further flatten, forming the modern alluvial plain and continental shelf basin.

Affected by the new Cathaysian tectonic system, the bedrock fault structure in the area was developed under a high activity intensity, and the main directions of the fault development were NNE, NE, and EW, followed by NW. The fault structures are all hidden type, and the major deep heat controlling faults include the Qihe-Guangrao fault, the Liaocheng-Lankao fault, and the Cangdong fault [28].

2.2. Stratigraphy

The North Shandong Plain is a vast and flat plain, covered by Quaternary sediments, with bedrock outcrops in some local areas, such as the Wudi Mountain. Since the Miocene, the plain has been slowly sinking, and thick Cenozoic strata have been deposited, including the Paleogene Kongdian, Shahejie, and Dongying Formations; the Neogene Guantao and Minghuazhen Formations; and the Quaternary Pingyuan Formation. The upper part of the Guantao Formation is composed of fine-medium sandstone and mudstone interbedded with siltstone, while the lower part is composed of thick bedded massive conglomerate and pebbled sandstone. Its bottom layer was widely developed, with sandy conglomerate containing quartz and flint.

2.3. Geothermal Reservoirs

In the North Shandong Plain, the geothermal water resources are mainly enriched in the Paleogene and Neogene porous bedded sandstone. The main geothermal reservoirs include the lower section of the porous Minghuazhen Formation and the porous Neogene Guantao Formation and the porous Paleogene Dongying Formation. Among them, the Guantao Formation is the most valuable reservoir and has the largest geothermal potential (Figure 3). The reservoir occurs throughout the plain, except in the eastern part of the Guangrao County [5, 29]. Its burial depth gradually deepens from the southwest to the northeast, generally exceeding 1,000 m in the depression, with a thickness ranging from 200 m to 700 m. The main reservoir is the coarse-grained sandstone and glutenite in the lower parts of the Guantao Formation (Ng2 in Figure 3). The geothermal well temperature measurement data show that the temperature of the water in the Guantao Formation is generally 50°C–90°C. It is a heat conduction dominated geothermal resource.

3.1. Study Methods

Based on temperature logging, water level measurement, well testing, hydrochemical analysis, and radioactive isotope analysis data for 800 geothermal boreholes, in this study, the geotemperature field, paleosedimentary facies, geothermal water flow field, hydrochemical field, and geothermal water age and abundance were analyzed. Then, the sources of heat and water and their enrichment were determined.

In general, the thermal conductivity of dense hard rock is high, while those of the soft silty sand and soft soil are low. For example, the thermal conductivity of limestone and granite is 2.01–2.721 W/(m·°C), that of sandstone is ~2.596 W/(m·°C), that of loose sand is 0.264–0.712 W/(m·°C), and that of sandy clay is ~0.921 W/(m·°C). Groundwater activity is mainly reflected by the migration speed. Generally, the faster the migration speed, the lower the geothermal gradient, and the slower the migration speed, the higher the geothermal gradient. Therefore, based on vertical temperature measurement data and borehole data from geothermal wells, the correlations between the lithology, groundwater activity, and geothermal gradient can be explored [3033].

The hydrochemical characteristics of geothermal water include the hydrochemical types, component concentrations, and isotopic features. The hydrochemical type of geothermal water is related to the main ions in the geothermal water. Based on the hydrochemical types of water samples and their distributions, as well as the geological characteristics of water samples, a certain rule can be obtained. Isotopic analysis is a key method of investigating the geothermal water formation mechanism and its recharge, runoff, and drainage conditions [3436]. Methodologies based on the use of isotopes in a full spectrum of hydrological problems encountered in water resource assessment, development, and management activities are already well established and are an integral part of many water resource investigations and environmental studies [37].

The isotopic analysis mainly includes the 14C and hydrogen and oxygen isotope characteristics. The hydrogen and oxygen isotope characteristics are the most widely used isotopic analysis method in geothermal water isotope research. By comparing the specific values of the hydrogen and oxygen isotopes with the China atmospheric precipitation line, most researchers currently believe that atmospheric precipitation is the recharge source of geothermal water. 14C isotope analysis is mainly used to determine the age of ancient geothermal water. The flow direction and evolutionary migration path of geothermal water can also be determined based on the age difference, the regional geological structure, the paleosedimentary facies, and the mineralization degree of the geothermal water [3841].

The analysis of the enrichment mechanism of geothermal water includes the characteristics of the geothermal reservoir’s paleosedimentary facies and geothermal water flow field and its enrichment rules. The basal conglomerate is the main mark and the material basis of the water in the Guantao Formation reservoir. Based on the thickness, burial depth, and distribution features of the basal conglomerate of the Guantao Formation obtained from the borehole data, as well as the regional geological structure and hydrogeological characteristics, the paleosedimentary facies of the Guantao Formation geothermal reservoir can be identified, and the distribution and enrichment of the geothermal water resources can be understood, which is of great significance for studying its formation mechanism and for scientific development of the Guantao Formation geothermal reservoir [4246].

According to the geothermal water level contour map and the burial depth zonation map drawn from the water level measurement data for the geothermal wells in different periods, the characteristics of the geothermal water flow field, the migration direction, and the water level fluctuation amplitude in different periods and different regions can be analyzed. By projecting the single geothermal well yields onto the map and drawing a zonation map, the water enrichment rules can be clearly determined. In addition, the relationship between the paleosedimentary facies and the water abundance can be determined.

3.2. Sampling and Analysis

Fifty-three geothermal water samples for isotopic analysis and 13 samples for chemical analysis were collected from geothermal wells between 2016 and 2017. Their locations and data are presented Tables 1 and 2. Prior to sampling, water was pumped for more than 30 min in order to obtain fresh water from the reservoir aquifer. The pH was measured in the field using a portable pH meter (HACH HQ40d) that had been calibrated beforehand. The alkalinity was determined via titration with 0.05 M HCl on the day of the sampling. The fresh water samples pumped from the reservoir aquifer were filtered on site through 0.45 μm membranes prior to collection in prewashed polyethylene bottles.

The samples for cation analysis were acidified (with ultrapure HNO3 to pH=2). The anions were determined using ion chromatography (Dionex ICS-1100), and the cations were determined using Inductively Coupled Plasma Optical Emission Spectroscopy (ICP-OES) (Thermo Fisher ICAP-6300) at the Institute of Hydrogeology and Engineering Geology of Shandong Province within 2 weeks of sampling. The charge balance errors of the water samples were generally less than 5%, which is within the limits of acceptability.

The stable isotope compositions (2H and 18O) of the water samples were measured using the LGR IWA-35-EP at the Institute of Hydrogeology and Environmental Geology, Chinese Academy of Geological Sciences. The δ2H and δ18O values are reported relative to the Vienna Standard Mean Ocean Water (VSMOW) in the conventional δ (‰) notation, with precisions of ±0.6 and ±0.2‰, respectively.

The 14C activities of the dissolved inorganic carbon (DIC) were measured radiometrically using a liquid scintillation spectrometer (PE 1220 QUANTULUS) after conversion to benzene at the Institute of Hydrogeology and Environmental Geology, Chinese Academy of Geological Sciences, to determine the residence time of the geothermal groundwater. The 14C results are reported in standard notation as percent modern carbon (pMC). The radiocarbon ages were calculated using a 5730-year half-life for 14C and are reported in years before present (BP), with 1950 taken as the starting year. The accuracy of the age estimates could possibly be improved by using a correction scheme for the matrix carbon [47]. The δ13C ratio of the DIC was subsequently analyzed using an isotope ratio mass spectrometer and liquid scintillation counting after conversion to benzene. The δ13C ratios are reported relative to Pee Dee Belemnite (PDB), and the 14C is reported using the unit of pMC.

3.3. Data Analysis

In this study, temperature measurement, water level, well yield, and water quality data were collected from 301 boreholes with depths of 1300–1800 m in 38 counties in the North Shandong Plain from 1998 to 2017. The temperature measurements included downhole and wellhead measurements. The downhole data were measured using a SKD-3000B type logging instrument when the well was completed. The wellhead data were measured during the pumping tests. The geothermal water level included the monitoring data collected in 2005 and 2017. The water level in 2005 represents the initial natural state before large-scale development, while that in 2017 represents the current situation of heavy exploitation.

A number of bore holes used for the lithology, pumping tests, temperature logging, water level monitoring, hydrochemistry, and isotope data collection in the North Shandong Plain are listed in Table 3, and the geothermal water isotope analysis data are presented in Table 1.

Based on the analysis of the data reliability, the statistical cluster analysis method was used to categorize all of the data and to average the sample values in the same area. Based on the regional distribution of and variations in the existing sample data and the geological conditions, in this study, the heat sources, water sources, and enrichment mechanism of the porous sandstone geothermal reservoirs were determined.

4.1. Geotemperature Field Characteristics

The formation of geothermal water in a region is affected by the mantle structure, the magmatic activity in the deep part of the crust, the geological structure, the stratigraphic lithology, and the groundwater activity in the shallow part of the crust. In addition, it is related to factors such as the caprocks (thermal insulation layer), geothermal reservoir space, heat sources, and geothermal water supply sources [28, 4851].

The caprock of the Guantao Formation geothermal reservoir is a soft layer composed of clay and sand in the loose sedimentary layers of the Quaternary Pingyuan and Minghuazhen Formations. The sandstone and sandy conglomerates of the Guantao Formation are rich in geothermal water, which not only transfers heat but also stores heat and water. This lithology provides a good storage space, constituting a porous stratiform geothermal reservoir. The heat sources are mainly heat conducted from the normal deep crust and upper mantle. The large and deep ultracrustal faults, such as the Cangdong, Liaocheng-Lankao, and Qihe-Guangrao faults, play an important role in the communication with and transfer from the magma heat source in the deep crust and upper mantle, and they constitute good channels for deep heat flow upwelling.

4.1.1. Caprock Geothermal Gradient Characteristics

The crustal surface geotemperature characteristics, which are a reflection of the thermal energy change in the interior of the Earth, are controlled by the regional geological development history and the geological structures. The main factors that influence the geotemperature include the terrestrial heat flow, basement fluctuations, stratigraphic lithology, fracture activity, magmatic activity, and groundwater movement. Based on the temperature measurement data, the geothermal gradient is 3.0–4.0°C/100 m (Figure 4). The horizontal distribution characteristics of the geothermal gradient are as follows.

  • (i)

    The distribution of the caprock’s geothermal gradient is clearly controlled by the distribution of the geological structures. In the positive tectonic zone, the burial depth of the bedrock is shallow, and the geothermal gradient is large, ranging from 3.6 to 4.0°C/100 m. In the negative tectonic zone, the burial depth of the bedrock is deep, and the geothermal gradient is small, generally less than 3.4°C/100 m. In the uplift area of the shallow buried Archaeozoic crystalline basement and the Mesozoic volcanic rocks, the geothermal gradient is generally greater than 3.6°C/100 m. The geothermal gradient contour line is basically consistent with the NNE-trending basement structure in the study area

  • (ii)

    The Mesozoic igneous rocks, especially the intrusive rocks, significantly influence the geothermal gradient. In the area with shallow buried Mesozoic igneous intrusive rocks, such as in the northwestern corner of Jinan in the Zouping-Huantai region, the geothermal gradient is generally greater than 3.6°C/100 m

  • (iii)

    The geothermal gradient in the depression area adjacent to the basement uplift area is quite low, such as in the Yangxin-Wudi area. Under the influence of the Chengning uplift in the western part of the study area, the uniform heat flow from the deep part of the Earth is concentrated toward the Chengning uplift area, resulting in a decrease in the geothermal heat flow in this area and a corresponding decrease in the caprock’s geothermal gradient

  • (iv)

    The deep and large Liaocheng-Lankao and Qihe-Guangrao faults (Figure 2), which cut into the upper mantle, significantly control the geotemperature field. The high geothermal gradient contours are distributed along the strike of the fracture belt

4.1.2. Temperature Distribution Characteristics of the Geothermal Reservoir

A zonation diagram of the temperature of the geothermal reservoir at a depth of 1,500 m in the different regions was drawn based on the single geothermal well temperatures and their depths (Figure 5).

As shown in Figure 5, the temperature at a depth of 1,500 m in the North Shandong Plain ranges from 55°C to 70°C. The high-temperature zones are mainly located in the uplift zones, such as the Guanxian, Gaotang, Gaoqing-Boxing, Ningjin-Qingyun, Zhanhua-Kenley, and Gudao uplifts. Although some local areas are not uplift regions, they are at the intersections of faults, where the tectonic structures communicate with the deep heat source. It can be seen that the geotemperature of the Guantao Formation in the North Shandong Plain is controlled by the structural framework, with high temperatures at the intersections of faults and in uplift areas and relatively low temperatures in other areas.

4.1.3. Vertical Variations in the Geothermal Gradient

The spatial distribution characteristics of the regional geotemperature field can be expressed in the form of the geothermal gradient. Based on the geothermal well drilling and temperature logging data, the vertical variations in the geothermal temperature are as follows. The geothermal temperature increases with increasing burial depth in the constant temperature zone. The vertical variation in the geothermal gradient is mainly affected by the thermal conductivity of the rock and the groundwater movement. When the conductivity is low, the geothermal gradient is high; otherwise, the geothermal gradient is low. In addition, when the groundwater activity is strong, the geotemperature gradient is low.

The general features are as follows. The Quaternary layer, which has a loose structure, has a low conductivity, a large thermal resistance, and a high geothermal gradient; however, it is influenced by the groundwater activity; so, the geothermal gradient is reduced. The Paleogene and Neogene layers, which has a denser structure, have a higher conductivity than that of the Quaternary layer due to the basically static state of the groundwater stored in the formation and the small influence of the groundwater convection on the heat energy; this layer has a higher geothermal gradient. The bedrock, which has the densest structure, has a higher conductivity and a correspondingly low geothermal temperature gradient (Figure 6).

4.2. Hydrochemistry Characteristics of Geothermal Water

4.2.1. Hydrochemistry Type of Geothermal Water

According to the analysis of the water quality data from the geothermal wells (Table 2), the main hydrochemistry Schukarev type of the Guantao Formation geothermal water is Cl-Na, followed by Cl·SO4-Na and SO4·Cl-Na. The Cl-Na type geothermal water is mainly distributed in the Linqing and Dezhou-Linyi depressions and their vast eastern area. The Cl·SO4-Na type geothermal water is mainly distributed in the Guanxian-Shenxian-Wucheng-Pingyuan area, and the SO4·Cl-Na type is mainly distributed in eastern Liaocheng-Gaotang. The horizontal variation characteristics of the hydrochemistry types are as follows: the hydrochemistry type gradually transitions from SO4·Cl-Na type in the southwest to Cl·SO4-Na type in the northwest and north and to Cl-Na in the east and northeast (Figure 7). It is speculated that desulfurization could occur during this process.

4.2.2. Major Constituents of Geothermal Water and Their Characteristics

As shown in Table 2, the main cation in the Guantao Formation geothermal water is Na+, with concentrations of 1,132.5–6,187.37 mg/l, accounting for 74.95–91.7% of the total cations. The main anion is Cl-, with concentrations of 858.94–10,900.00 mg/l, accounting for 37.03–98.49% of the total anions. The TDS values are 4.07–18.52 g/l, and the pH values are ~7.14–8.1. Based on the TDS and pH values, the water is weakly alkaline. The TDS value gradually increases from west to east. Dongying, the central city located in the Yellow River delta, has the highest TDS value. The Na/Cl value is in the range of 0.81–2.15 with a mean of 1.25, which means a weak cationic exchange. The Na/K ratios are 36.45–1,857.14, with an average of 130, which are higher than that of normal groundwater. The Mg/Ca ratios are 0.03–0.79, generally about 0.3, which are lower than that of ground water with a normal temperature. Therefore, the geothermal water is depleted in potassium and magnesium. This may be the result of the transition from montmorillonite to chlorite.

4.2.3. Hydrochemistry Evolution Path of Geothermal Water

Cl- is a conservative ion, and it is generally difficult to precipitate it from groundwater [5254]. Therefore, the evolution path of geothermal water can be preliminarily inferred based on the variations in the Cl- content. As shown in Figure 6, in the western part of the North Shandong Plain, the Cl- concentration increases gradually from the western uplift zone to the Linqing and Jiyang depressions. Thus, it can be inferred from this phenomenon that the geothermal water may have undergone a certain degree of hydro-rock reaction or evaporation concentration from the western uplift zone to the Linqing and Jiyang depressions. The Cl- concentration of the geothermal water in the Chengning uplift gradually increases from north to south, and its concentration is higher than that on its both sides in the eastern region (Laoling—north of Zhanhua) and the western region (Linqing depression). Similarly, it can be speculated that the geothermal water in this area received peripheral lateral runoff recharge in the northern part of the Chengning uplift. The Cl- and Na+ concentrations of the geothermal water are higher than those in the surrounding area at the center of every depression in the North Shandong Plain, such as the Shenxian, Linqing, Linyi, Zhanhua, Huimin, and Dongying depressions. This phenomenon reveals that the evolution of the geothermal water’s migration path is from uplift to depression, indicating that the syngenetic water during the sedimentation is one of the sources of the geothermal water in the Guantao Formation, which later experienced intensive evaporation and concentration.

4.2.4. Isotopes and Their Characteristics

Based on the hydrogen and oxygen isotope data for the geothermal water in the Guantao Formation in the North Shandong Plain (Table 1) and the location of these data with respect to China’s meteoric water line (CMWL) (δD=7.7δ18O+7.5), the genesis of the geothermal water in the Guantao Formation in the North Shandong Plain can be determined (Figure 8).

Figure 8 shows that the D and 18O values of the geothermal water in the Guantao Formation basically plot close to the CMWL. This indicates that the source of the geothermal water is atmospheric precipitation. The δD and δ18O values of the geothermal water plot slightly below China’s atmospheric precipitation line, which indicates the dilution effect of steam water from the lower strata during the long geological history. During the runoff of the geothermal water from southwest to northeast, interactions occur between the geothermal water and the surrounding rock, and thus, the δ18O value increases, and obvious oxygen isotope drift occurs.

According to the corrected 14C data (Table 1), the age of the geothermal water in the North Shandong Plain is generally 5–25 kyr. The age of the geothermal water in the depression area is older, ranging from 15 to 25 kyr, while the geothermal water in the uplift region is younger, ranging from 5 to 20 kyr (Figure 9). This confirms that the evolution of the geothermal water’s migration path is from the uplifts to the depressions.

4.2.5. Water–Rock Interactions

The Na-K-Mg ternary diagram was proposed by Giggenbach in 1988 and is commonly used to evaluate the water–rock equilibrium state and to distinguish between different types of water samples [55]. It is a graphical method that can be used to distinguish geothermal water from completely balanced water, partially balanced water, and immature water. It can be used to determine the types of geothermal water samples, evaluate the equilibrium state, and forecast mixing trends between hot water and the surrounding rock. Figure 10 shows the Na-K-Mg ternary diagram for the geothermal water samples from the Guantao and Dongying Formations in the North Shandong Plain. As can be seen from Figure 10, the hot water in the sandstone of the Guantao Formation and the Dongying Formation is in the local equilibrium area, i.e., immature water. It plots very close to the Mg end member, indicating that the water–rock interactions are weak, the ionic equilibrium between the water and rocks has not been reached, and dissolution is still occurring.

4.3. Enrichment Mechanism of Geothermal Water

Based on the geothermal water’s occurrence environment, the characteristics of the geothermal water flow field, and the water abundance, the rule of the geothermal water enrichment was investigated.

4.3.1. Lithology of the Geothermal Reservoir

According to the data from 301 geothermal wells, the main lithology of the geothermal reservoir is the fine sandstone, coarse sandstone, conglomeratic sandstone, and glutenite in the fluvial and alluvial fan facies. The conglomerate grains are semicircle shaped, with medium roundness. The sandstone/glutenite reservoir accounts for 30–40% of the thickness of the Guantao Formation, with a single layer thickness of 2–40 m (Figure 11).

The lithology of the geothermal reservoir is characterized by a fining upward positive cycle (i.e., fine particles in the upper part and coarse particles in the lower part), with basal conglomerate. The composition of the gravel is dominated by quartz and black chert, with a diameter of 1–10 mm and medium roundness. The conglomerate is characterized by weak cementation and poor diagenesis, resulting in an unconsolidated structure. The loose structure and coarse grains lead to a high porosity of 20–40% and a high permeability of 100–1500 mD (Figure 12). Thus, the Guantao sandstone/glutenite reservoir has good water storage space. In the horizontal direction, the lithology is characterized by coarse grains in the central-western and eastern regions and fine grains in the southern and central-eastern regions.

4.3.2. Porosity versus Permeability

The effect of grain size on the permeability has been determined for unconsolidated sand [56] and for sandstones. At the time of deposition, medium-grained sandstones have a higher permeability than very fine- and fine-grained sandstones [57, 58]. The average permeability also gradually decreases with a decreasing degree of sorting for both unconsolidated sand [56] and sandstones [57, 59]. Several factors affect the porosity and permeability of the sandstone/glutenite of the Guantao Formation in the North Shandong Plain. In addition to burial depth, the grain size also has a strong influence on the permeability [14]. The burial depth of the Guantao Formation in North Shandong Plain is generally around 1000 m. Thus, the grain size is the most important factor influencing the porosity and permeability in this area. As shown in Figure 12, the relationship between the porosity and permeability changes systematically with grain size; so, for a given porosity, fine-grained and medium-grained sandstones have lower permeabilities than coarse-grained sandstones and glutenite; that is, at a given porosity, the coarse and thick sandstone/glutenite reservoirs in the ancient river channels have higher permeabilities than the fine and medium grained reservoirs. This relationship is consistent with those of other sandstone reservoirs around the world. For example, Weibel et al. concluded that the shallow buried (<2500 m) fluvial and upper shore facies sandstones of the Gassum Formation are typically good geothermal reservoirs in the Danish area [14]. The permeabilities of these shallow buried sandstones vary with grain size, and for a given porosity, the permeability increases with increasing grain sizes (from very fine- to fine- to medium- to coarse-grained).

4.3.3. Geothermal Water Occurrence Environment—Paleosedimentary Facies of the Guantao Formation Geothermal Reservoir

The most important feature of the ancient river channel is its associated riverbed facies, which is dominated by pebbles or coarse sand at the bottom and transitions to sand or silty sand at the top. In a vertical section, the particle size changes from coarse at the bottom to fine at the top. In a longitudinal section, the upstream area is thicker than the downstream area. Buried paleochannels are generally rich in groundwater. Studying paleochannels helps us to understand the historical evolution of rivers and landforms in the region, and it is an important means of studying the ancient geographical environment.

The formation age of the geothermal water in the Guantao Formation reservoir in the North Shandong Plain is relatively old. The reservoir is classified as a layered porous sandstone geothermal reservoir. The basal conglomerate is the most obvious indicator of the fluvial sedimentary environment in the paleochannel of the Guantao Formation, and it is also the material basis for the abundance of geothermal water in the Guantao Formation reservoir. Therefore, based on the large amount of data on the distribution of the basal conglomerate and the burial depth and thickness revealed by the geothermal boreholes (Table 4), we identified the paleosedimentary facies of the Guantao Formation in the North Shandong Plain to delineate the distribution range of the paleochannel zone (Figure 13) and to determine the geothermal water enrichment zone.

As can be seen from Figure 13, the burial depth of the Guantao Formation in the North Shandong Plain is generally greater than 1000 m, except for in the local uplift areas, and the thickness of the basal conglomerate is generally 50–150 m. The sedimentary facies include two large paleochannel zones. (i) The Linqing-Wucheng-Dezhou-Lingxian paleochannel zone, which is controlled by the Gaotang uplift in the south and the Chengzikou-Ningjin uplift in the east, forms a large alluvial-pluvial fan with coarse and thick gravel stone deposits. (ii) The Zhanhua-Hekou-Xianhe paleochannel zone, which is controlled by the Diaokou and Yihezhunang uplifts in the northwest and the Wudi uplift in the south, also forms a large alluvial-pluvial fan with coarse and thick gravel stone deposits. Influenced by the Chenzhuang and Qingtuo uplifts in the downstream area, the paleochannel was bifurcated into two belts, with the southern part distributing along the Lijin-Dongying.

4.3.4. Characteristics of Geothermal Water Flow Field

In the early part of this century, the exploitation degree of the geothermal water in the Guantao Formation was very small; so, the geothermal water was in a natural self-flowing state, and the flow field was basically stable. According to the water level data reported during the completion of the geothermal wells earlier in this century (Table 5), the general flow direction of the geothermal water in the Guantao Formation was from west to east or from southwest to northeast.

In the last 10 years, as the number of new geothermal wells and the mining scale have increased year by year in the urban areas in various counties, geothermal water depression cones have been formed. These depression cones are centered in the urban areas, such as in Liaocheng, Dezhou, Linqing, and Pingyuan, and have changed the natural flow direction.

Although the overall runoff direction is still from southwest to northeast, the geothermal water migrates from the periphery to the center of the depression, and the hydraulic gradient also changes significantly.

The burial depth of the water level in the depression cones in some counties is generally greater than 30 m and that in urban Dezhou exceeds 70 m. The burial depth of the water level in the noncone areas is ~20 m. The maximum decrease in the geothermal water level over the years is more than 80 m (Figures 14 and 15).

4.3.5. Geothermal Water Enrichment

According to the pumping test data for the geothermal wells in the Guantao Formation collected when the wells were completed, the single well yield was generally in the range of 60–100 m3/h (Table 6), and the water-rich areas were mainly distributed in the Linqing-Wucheng-Dezhou belt in the northwest and in the Dongying-Hekou-Xianhe belt in the east, with single well yields of greater than 85 m3/h. The areas with poor water abundance conditions were mainly distributed in the Liaocheng and Binzhou-Boxing region in the southern and central parts of the North Shandong Plain, with single well yields of less than 70 m3/h, while the yields in other areas were generally in the range of 70–85 m3/h (Figure 16).

As shown in Figure 16, the regularity of the geothermal water enrichment in the Guantao Formation is highly consistent with the distribution area of the ancient river channel belt; that is, the water-rich area is distributed in the ancient river channel area. The closer the water-rich area is to the ancient river channel belt, the deeper and coarser the buried conglomerate aquifer is, and the more water-rich it is, and vice versa. The most water-rich conditions occur in the alluvial fan in the ancient channel belt.

  • (i)

    The reservoir temperature in the uplift areas and its intersection with the faults is high, while that in other areas is relatively low. The highest temperatures mainly occur at the intersection of the uplifts and the faults. Vertically, the geothermal gradient is comprehensively influenced by the thermal conductivity of the rocks and groundwater activity. The lower the thermal conductivity, the higher the geothermal gradient

  • (ii)

    The chemical types of the geothermal water in the study area are mainly Cl-Na, followed by Cl·SO4-Na and SO4·Cl-Na. The TDS values of the geothermal water are 4.07–18.52 g/L, and the pH values are 7.14–8.1. In general, the TDS and Cl- contents gradually increase from west to east and from the local uplifts to the depressions, which also reflect the flow path of the geothermal water

  • (iii)

    The relationship between the porosity and permeability changes systematically with grain size. For a given porosity, the fine-grained and medium-grained sandstones have lower permeabilities than the coarse-grained sandstones and glutenite

  • (iv)

    The single well yield is generally 60–100 m3/h, and the distribution of the geothermal water enrichment is consistent with that of the ancient river channel belt; that is, the water-rich area is mainly distributed in the two ancient river belts, and the single well yields in these areas are generally greater than 85 m3/h. The closer to the ancient channel belt, the deeper and coarser the buried sand conglomerate aquifer is, and the more water-rich the aquifer is, and vice versa. The most water-rich area is the alluvial-proluvial fan in the ancient river channel belt

  • (v)

    The geothermal water is sourced from atmospheric precipitation. The corrected 14C ages are generally ~5–25 kyr, with the older geothermal water occurring in the depression area. The evolution of the migration path of the geothermal water is from the uplift area to the depression area, which is consistent with the evolutions of the paleosedimentary environment of the geothermal reservoir, the geothermal flow field, and the water quality

The borehole data and geochemistry data used in this paper are available from the authors upon request.

The authors declare that there is no conflict of interest regarding the publication of this article.

This work was supported by the National Natural Science Foundation of China (grant numbers U1906209 and 41877192).

Exclusive Licensee GeoScienceWorld. Distributed under a Creative Commons Attribution License (CC BY 4.0).