Characteristics of Deep Groundwater Flow and Temperature in the Tertiary Pohang Area, South Korea

Featured Application

This study can reveal the hydraulic characteristics of deep groundwater in sedimentary and volcanic rocks.

Abstract

The geothermal characteristics of specific areas on Earth can be identified using geophysical and lithological logs based on deep boreholes, such as those more than 1000 m in depth. Based on the combined analyses of geophysical logs (temperature, caliper, electrical conductivity (EC), and natural gamma-ray logs) and lithological logs, as well as 1D steady-state heat transfer model, the deep groundwater flow and temperature were characterized in four deep boreholes (BH-1–BH-4) roughly 2000 m deep in the area of Heunghae-eup, Pohang, South Korea. The estimated thermal gradients from the temperature profiles are as follows: 22.37–30.77 °C/km for BH-1, 35.67–64.52 °C/km for BH-2, 40.85–46.44 °C/km for BH-3, and 33.33–35.71 °C/km for BH-4. According to the geophysical logs and lithology profiles, the groundwater mainly flows into and out of the boreholes through the basic dyke, rhyolite, and sandstone/mudstone. Evidently, the groundwater flows moving through the fractures and faults induce nonlinear temperature changes. The upward and downward groundwater flows passing through fractures and faults can be estimated using a 1D steady-state heat-transfer equation, by considering a fracture angle based on the lithological and geophysical profiles. To determine the direction (up/down) and rate of groundwater flow, the values for the parameter β were estimated as follows: −1.95 to 5.40 for BH-1, −13.48 to 4.87 for BH-2, −1.76 for BH-3, and −3.39 to 14.15 for BH-4.

1. Introduction

The geothermal flow of the Earth’s interior is generated mainly by conduction and convection. Based on conduction and convection, the geothermal temperature increases steadily with depth. The equilibrium temperature of groundwater in a borehole with an overall linear or nonlinear thermal gradient represents the Earth’s temperature at that depth. Heat transfer inside the Earth occurs predominantly through conduction. The temperature distribution in the continental crust and lithosphere is governed largely by heat generated from the radioactive decay of rocks and by heat conduction from mantle to the surface through heat-transfer media and pathways [1]. Thermal gradients have been estimated by using constant temperature gradients in states unaffected by fluid flows or based on the average slope of an entire section of boreholes from the surface and bottom temperatures (using temperature-log data) [2,3].

In South Korea, due to the lack of thermal energy generated from young volcanic activity [4], research on deep geothermal energy began only in the 2000s. The Pohang area in South Korea has the largest geothermal energy potential, with heat flows greater than 80 mWm⁻² [4,5]. Moreover, the Pohang area, with the youngest tectonic movements in the Korean peninsula, has possible deep-stressed fracture zones [6,7,8]. To evaluate the geothermal energy potential, the Korea Institute of Geoscience and Mineral Resources drilled four deep boreholes (BH-1 to BH-4) from January 2003 to November 2006. BH-4 is a pumping borehole for developing deep groundwater; it has the largest diameter (200 mm at the final depth) of the four boreholes, with a depth greater than 2 km [9]. The distance between BH-1 and BH-2 is 157 m, and the distance between BH-3 and BH-4 is 123 m. As geothermal water in Korea has a higher probability of geothermal flow along secondary pores, such as faults and fractures, than along primary pores, it is important to understand how deep geothermal water is associated with faults and geological structures [10]. The construction of a geothermal binary power plant using an enhanced geothermal system with a planned generation capacity of 1.5 MW was expected to be completed in the Pohang area by 2017. However, its construction was stopped owing to the Mw 5.4 Pohang earthquake on November 15, 2017, which damaged buildings in Pohang.

Bredehoeft and Papadopulos [11] derived a simple steady-state solution for estimating vertical water fluxes from deep temperature profiles in homogeneous environments. This approach has been used in several studies, i.e., to calculate water fluxes from the curvatures of temperature profiles recorded in boreholes [12,13]. Furthermore, the estimation of the borehole vertical flow is of specific interest, as it can be used to derive large-scale hydraulic connections [14]. In addition, the temperature profiles display a strong correlation with the vertical borehole velocity profiles [15].

The aim of this study was to characterize deep groundwater flows and temperatures, based on the geophysical logs and lithology profiles of four deep boreholes in Heunghae-eup, Pohang, South Korea (Figure 1). The lithology profiles and geophysical logs can be used to effectively detect the fractures intersecting a borehole. Additionally, the temperature, caliper, electrical conductivity (EC), and natural gamma ray logging can be implemented for geological and stratigraphic interpretation, thermal gradient measurement, groundwater flow pattern determination, casing installation, and screen section determination [16,17].

2. Geological Setting

Pohang belongs to the East Sea alkaline volcanic zone [18], occupying the central part of the Tertiary Pohang Basin. Granodiorite bedrock was formed as a result of the igneous activity at the end of the Paleozoic Era and is overlain by repeated layers of Cretaceous sandstone/mudstone (Figure 1 [19]). Cretaceous lapilli tuff overlies the Cretaceous sedimentary rocks. Layers of Eocene crystal tuff and lithic crystal tuff are present on the Cretaceous rocks. In the middle Eocene, biotite granite intruded on the older rocks. Moreover, Miocene semiconsolidated sedimentary rocks (the Cheonbuk conglomerate, Hakrim Formation, Heunghae Formation, Idong Formation, Duho Formation, and Yeonam Formation) belonging to the Yeonil Group and composed of sandstone, shale, and mudstone were deposited along the east/southeast direction as a result of several marine transgressions and regressions, with a thickness greater than 1000 m [20,21,22,23].

The four deep boreholes were drilled vertically in the study area, at 1.8–2.8 km from the epicenter of the Mw 5.4 Pohang earthquake in 2017 (Figure 1). The borehole diameters decrease with depth. BH-1 has a depth of 1100 m, with a diameter of 76.2–101.6 mm; BH-2 has a depth of 1505 m, with a diameter of 88.9–203.2 mm; BH-3 has a depth of 950 m, with a diameter of 254–355.6 mm; and BH-4 has a depth of 2383 m, with a diameter of 203.2–254 mm. From top to bottom, the rock types of BH-1, BH-2, BH-3, and BH-4 consist of semiconsolidated mudstone, crystal tuff, lapilli tuff, sandstone/mudstone, rhyolite, andesitic breccia, dacite, and granodiorite; these were identified using a drilling slurry. The thickness of the rock types is as follows: 200–370 m for the semiconsolidated mudstone, 73–240 m for the crystal and lapilli tuffs, 46–907 m for the mudstone/sandstone, 259–375 m for the rhyolite, and 834 m for the andesitic breccia [9,24].

The depth to the water is −60 m at BH-1 (−45 m, mean above sea level (masl)), −120 m at BH-2 (−107 m, masl), −120 m at BH-3 (−109 m, masl), and −110 m at BH-4 (−99 m, masl).

3. Methods

3.1. Theoretical Background of Thermal Heat Transfer

The geothermal flux inside the Earth variably affects the thermal gradient [25]. The transient heat transfer in homogeneous and isotropic saturated porous media can be expressed as follows [26]:

$$\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}+\frac{{\partial }^{2}T}{\partial z^2}-\frac{c_0{\rho }_0}{\kappa }\left[\frac{\partial \left(v_{x}T\right)}{\partial x}+\frac{\partial \left(v_{y}T\right)}{\partial y}+\frac{\partial \left(v_{z}T\right)}{\partial z}\right]=\frac{c\rho }{\kappa }\frac{\partial T}{\partial t}$$

(1)

where T denotes the temperature at time t; $c_0$ and ${\rho }_0$ denote the specific heat and density of the fluid, respectively; c, ρ, and κ denote the specific heat, density, and thermal conductivity of the rock matrix in fluid-filled rocks, respectively; and v_x, v_y, and v_z denote the fluid velocities in the x-, y-, and z-directions, respectively. In the case of 1D heat conduction in the vertical direction [11], the assumptions are as follows: the borehole is filled with fluid at a steady state and without the influence of drilling and pumping; thermal conductivity is constant at the object section; fluid flow occurs only in the vertical direction; temperature changes in the borehole occur owing to the vertical flow; and the fluid flow rate of the object section is uniform. From Equation (1), the 1D steady-state heat transfer equation for an inclined fracture zone can be expressed as follows [11]:

$$\frac{{\partial }^{2}T}{\partial z^2}-\frac{c_0{\rho }_0}{\kappa }\frac{\partial \left(vT\right)}{\partial z}=0$$

(2)

where v is the fluid velocity vector. Equation (1) can be simplified as Equation (4) based on Ge’s assumptions [27], as follows:

$$\frac{\partial \left(vT\right)}{\partial z}\approx v_{\xi }\gamma$$

(3)

$$\frac{{\mathrm{d}}^{2}T}{\mathrm{d}{z}^2}-\frac{\mathsf{\alpha }\gamma }{L}=0$$

(4)

where v_ξ denotes the fluid velocity vector in the fracture direction (the direction of the local coordinate ξ), γ denotes the average magnitude of the temperature gradient in the fracture zone, and $\alpha$ denotes the Peclet number, equal to $\frac{c_0{\rho }_0{v}_{\xi }L}{\kappa }$.

When the upper and lower boundary temperatures of the aquifer or geological formation are T = T₀ (at z = 0) and T = T_L (at z = L), respectively (Figure 2), Equation (4) can be expressed as follows [27]:

$$\frac{T_z-T_0}{T_L-T_0}=z_D+\frac{\beta }{2}{z}_D\left(z_D-1\right)$$

(5)

where L denotes the thickness of the fractured aquifer, i.e., the vertical distance between z = 0 and z = L. The parameter z_D (dimensionless) is equal to z/L. The parameter β (dimensionless) can be expressed as follows:

$$\beta =\mathsf{\alpha }\sin \omega$$

(6)

where $\omega$ denotes the average angle of the fracture zone from the horizontal plane. The parameter β has a different sign depending on the direction of the flux: zero in the case of conduction, positive in the case of downward flow, and negative in the case of upward flow.

3.2. Determining Thermal Conductivity, Specific Heat, and Density

Thermal conductivity is a vital factor for determining the temperature inside the Earth [28,29]. In this study, a Nanoflash instrument (model LFA-447 of Netzsch Co., Germany) equipped with an infrared thermometer was used to measure the thermal diffusivity, by transferring heat from the bottom side to the top side of the rock sample. The LFA-447 was calibrated with a standard sample (pyroceram 9606) of thermal diffusivity of 1.926 mm²/s at 25 °C; the measurement error is no greater than ±3% between the measured value and true value.

The thermal conductivity of a dried rock sample can be estimated as follows:

$$\kappa =\rho ·D·c_p$$

(7)

where ρ denotes the bulk density of the sample, D denotes the thermal diffusivity, and c_p denotes the specific heat. The dry sample is prepared in a vacuum dryer for 24 h.

The thermal conductivity changes significantly according to the dried/saturated state and porosity of the sample. The weighted geometric mean of the thermal conductivity of a saturated sample can be expressed as follows [30]:

κ_sat = κ_water ⁿ·κ ¹⁻ⁿ

(8)

where κ_water denotes the thermal conductivity of water (0.62 W/mK), κ denotes the thermal conductivity of the dry sample, and n denotes the porosity (used to determine the thermal conductivity of the saturated natural rock). The porosity was determined from the volume differences between saturated and dried samples. Yalcin and Dumlu [31] assumed a constant specific heat of 4.2 kJ/kgK and density of 1000 kg/m³ for estimating the vertical flow velocity. In contrast, in this study, specific heat and density values were applied for different temperature intervals of the temperature profile in the boreholes, as shown in

Table 1. Estimated vertical velocity, dip of fracture/fault, and thermal properties in the boreholes. Negative and positive signs of flow velocity indicate upward and downward flows, respectively.

Borehole	Gradient (°C/km)	Depth Interval (m)	Lithology	T0–TL (°C)	Specific Heat of Water (kJ/kgK)	Density of Water (kg/m3)	Thermal Conductivity (W/mK)	β	Flow Velocity (m/s)	Flow Direction	RSS	Dip of Fracture/ Fault (°)
BH-1	22.37–30.77	358–410	Dacite	22.42–24.32	4.180–4.181	997.2–997.7	2.73	−1.61	−2.27 × 10−5	Upward	0.037	62.81
		488–514	Sandstone/mudstone	27.06–28.39	4.179	996.2–996.5	2.97	−1.95	−5.93 × 10−5	Upward	0.012	64.28
		540–566	Rhyolite	30.00–31.30	4.178	995.3–995.7	3.87	−1.88	−7.51 × 10−5	Upward	0.009	63.38
		586–646	Rhyolite	29.96–34.69	4.178	994.2–995.7	3.87	5.40	1.09 × 10−4	Downward	0.113	50.42
		890–998	Sandstone/mudstone	48.25–52.77	4.181–4.182	986.8–988.8	2.97	−1.13	−7.86 × 10−6	Upward	0.009	72.89
BH-2	35.67–64.52	580–610	Rhyolite	42.68–43.77	4.179	990.8–991.2	3.87	4.87	3.16 × 10−4	Downward	0.285	28.64
		618–652	Rhyolite	43.78–45.57	4.179–4.180	990.0–990.7	3.87	−13.48	−5.16 × 10−4	Upward	1.713	45.90
		736–757	Basic dyke	45.18–45.60	4.18	990.0–990.2	2.86	−3.34	−1.41 × 10−4	Upward	0.125	51.02
		955–976	Sandstone/mudstone	51.58–53.93	4.182	986.2–987.3	2.97	2.05	2.39 × 10−4	Downward	0.188	17.12
		1178–1189	Sandstone/mudstone	57.45–59.41	4.184–4.185	983.5–984.5	2.97	−1.51	−5.90 × 10−4	Upward	0.011	9.67
		1257–1273	Sandstone/mudstone	60.88–62.12	4.186	982.1–982.7	2.97	3.90	4.69 × 10−4	Downward	0.032	24.67
		1287–1356	Sandstone/mudstone	62.76–65.06	4.186–4.187	980.5–981.8	2.97	−1.69	−1.97 × 10−5	Upward	0.162	64.43
BH-3	40.85–46.44	138–381	Mudstone, crystal tuff, and basic dyke	25.40–30.12	4.178–4.179	995.6–997.0	3.02	−1.76	−5.27 × 10−6	Upward	0.138	90.00
BH-4	35.71–48.78	303–350	Crystal tuff	26.98–28.63	4.179	996.1–996.5	3.00	14.15	2.17 × 10−4	Downward	2.013	90.00
		394–480	Crystal tuff	26.97–31.32	4.178–4.179	995.3–996.6	3.00	8.96	7.96 × 10−5	Downward	1.479	70.52
		567–600	Sandstone/mudstone	34.41–35.50	4.178	993.9–994.3	2.97	−3.39	−8.07 × 10−5	Upward	0.011	65.64
		873–967	Sandstone/mudstone	43.93–54.12	4.179–4.183	986.1–990.7	2.97	3.12	4.28 × 10−5	Downward	0.925	33.81

.

3.3. Geophysical Logs

Temperature logging is mainly used to estimate groundwater flows and the locations of aquifers and to identify different thermal gradients [32,33,34]. The thermal gradients can be classified into three categories: linear thermal gradients, wherein the fluid temperature in the borehole is in equilibrium with the underground temperature (Figure 3a); multiple thermal gradients, wherein the borehole consists of several strata with different thermal conductivities (Figure 3b); and nonlinear thermal gradients, wherein the fluid flows through fractures (Figure 3c,d) [35].

A caliper log uses a mechanical caliper, acoustic caliper, and sonar caliper to continuously measure a borehole diameter. Caliper log anomalies indicate fracture depths, based on the enlargement of the borehole diameter [36,37]. Inflows to and outflows from the borehole can be attributed to the bedding planes or fractures identified by the caliper log. The advantages of this approach are clearly demonstrated by the number of possible water-producing zones indicated by the caliper logs. In geothermal drilling activity, a caliper log is mainly used for detecting cavities in a well, measuring the wellbore volume to estimate the required cement volume, determining the placing of casing in a well, inspecting the casing for damage either during or after drilling throughout the lifetime of the well, and calibrating the well depth (as combined with the lithological and natural gamma-ray logs) [38,39].

EC logging detects the EC values of water in various aquifers, thereby providing a general indication of water quality and hydraulic information [40,41]. A systematic analysis of EC logging data can predict the location of possible hydraulically active zones in a borehole [42]. In general, the EC depends on the porosity, water temperature, and salinity, as it is controlled by ionic diffusion and dispersion, as well as by the moving water. Under the equilibrium condition, the ionic content at a given depth in a water-filled borehole effectively represents a balance of the various processes at that depth [43,44]. A depth with less saline water may represent the mixing point of formation waters, and high-salinity zones can be attributed to inflows of saline water into the borehole [45].

Natural gamma-ray logging is useful with respect to the stratigraphic identification of basic intrusive rocks, such as rhyolite. In addition, gamma-ray logging is an effective method for identifying felsic rock among the basalt formations in a borehole, using a linear relationship between the silica content in the volcanic rock and the natural gamma-ray intensity. Moreover, a natural gamma-ray log is often used to discriminate between Mesozoic sedimentary rock and Tertiary formations. For example, the presence of sandstone in semiconsolidated mudstone exhibits a lower gamma-ray intensity [46,47,48].

In this study, the geophysical logger of Robertson Geologging (UK) with the console of Micrologger II was used for temperature log with an accuracy ±0.5 °C, EC log with an accuracy ±2.5% at 500 μ S/cm, and caliper log with an accuracy ±6 mm.

4. Results

4.1. Lithological Profiles

Figure 4 shows the stratigraphic columns of the four deep boreholes; it is evident that crystal tuff lies below the semiconsolidated mudstone at the following depths: 206–428 m for BH-3 and 217–443 m for BH-4. Crystal tuff was not observed for BH-1 and BH-2. The crystal tuff exhibits a moderately weathered grade and light-green color in the upper section (212–270 m in depth), with oxidized pyrite in the disseminated state and a relatively fresh grade in the lower section, which is light gray [9].

Dacite can only be observed in BH-1 and is located between the semiconsolidated mudstone and crystalline/lapilli tuff at a depth of 359.1–401.4 m. It is dark gray and has a generally dense lithology. At a depth of 372–374 m, the dacite is moderately weathered; at 400.7–401.4 m, it is highly weathered.

Numerous basic dykes and tuffaceous sediments are distributed between the crystal tuff layers. The basic dykes intrude to a shallow depth, with sufficient horizontal elongation at the sandstone and lapilli tuff zones. The lapilli tuff, which is light green, exhibits a moderately weathered grade over the entire section. In BH-3, rhyolite is distributed at a depth of 630–900 m with a thin layer of sandstone/mudstone in the bottom section, and sandstone (525–630 m in depth) underlies the lapilli tuff in a different manner from that in BH-2, even though the two boreholes have a similar geological stratigraphy. The lapilli tuff in BH-2 lies between 352 and 425 m from the land surface with a thickness of approximately 357–427 m.

Sandstone/mudstone mostly alternate in the lower part of the lapilli tuff with a few intercalations. The sandstone is greenish-gray, and the mudstone is dark gray and reddish-brown; similar colors can be seen in BH-1, BH-2, and BH-3. Slightly weathered coarse-grained light gray tuffaceous sandstone can be observed in BH-4 at depths of 900–900.3 and 902.9–903 m. In BH-1, BH-2, BH-3, and BH-4, the sandstone/mudstone is dark gray in the upper parts and reddish-brown in the lower parts at a depth of 1000 m, suggesting an environmental change from oxidation to reduction. In BH-3, dark-gray mudstone is predominant, except for two reddish-brown zones in the lower part, at a depth of 428 m.

Dark greenish-gray andesitic breccia can be observed at a zone lower than that of the sandstone/mudstone beds in BH-4; it cannot be observed in BH-1, BH-2, and BH-3. Andesitic breccia generally contains dark-red volcanic blocks and is fresher than the tuff. Medium-grained granodiorite contains a large amount of biotite and a minor amount of hornblende; it is light/dark gray and generally has a darker color than that of biotite granite.

The detailed lithologies, fractures/faults, and casing installations of the four boreholes are presented in Appendix A.

4.2. Geophysical Log Profiles

The geophysical logging was conducted under steady-state conditions (without the influence of drilling and pumping) in October 2010 and December 2014, after the drilling of BH-1, BH-2, BH-3, and BH-4 from January 2003 to November 2006. Temperature, caliper, EC, and natural gamma-ray loggings were conducted to interpret the lithology, aquifer, and flow characteristics of the four boreholes (Figure 5, Figure 6, Figure 7 and Figure 8). Owing to the difficulty in acquiring highly reliable water temperature data immediately after drilling, the geophysical logging was executed at the stabilization of the groundwater temperature. In BH-1, only the temperature log data were available for the entire depth; the caliper and natural gamma-ray logs were measured up to 359 m, owing to a borehole collapse. The accuracy of the geophysical logging for BH-3 was uncertain with respect to the EC and temperature data, owing to the influence of drilling mud. In BH-4, geophysical logging was conducted to a depth of 2000 m, owing to temperature and pressure sensor limits [24].

An accurate temperature log is critical in determining the thermal gradient and thermal conductivity; however, this was difficult to achieve, owing to the high specific weight of the drilling mud used to maintain the borehole wall during excavation and because of the casing installed after the excavation. The temperature log data for BH-3 have low accuracy because the data were obtained during drilling and because of the effect of the drilling mud. The thermal gradients of the boreholes (as well as the groundwater flow rate and direction) were estimated from the temperature logs. The decreasing tendency of the temperature from the surface to a depth of 100 m indicates the influence of the surface air temperature on the underground temperature [49,50,51].

4.2.1. BH-1

Owing to the collapse of the borehole, only temperature logging was conducted for the entire depth; the caliper and natural gamma-ray loggings were performed only in the upper zone (Figure 5). The temperature log zones were identified at the following depths: 358–410 m for the dacite section, 488–514 and 890–998 m for the sandstone/mudstone section, and 540–566 and 586–646 m for the rhyolite section (Figure 5b,c). The highest temperature for BH-1 was 56.28 °C at a depth of 1092 m. Based on the temperature log, the thermal gradient of BH-1 was 30.77 °C/km at a depth of 258–420 m, 22.37 °C/km at a depth of 430–715 m, and 25.01 °C/km at a depth of 820–1092 m.

4.2.2. BH-2

Compared with the lithology data for the other boreholes, more caliper log anomalies were observed at the following depths: 576.4, 636.4, 704–773 m for rhyolite and basic dyke and 900, 960.75, 1065.55, 1087.25, 1124.8, 1181, 1126.9, and 1245 m for sandstone/mudstone (Figure 6a).

The temperature log zones were identified at depths of 580–610, 618–652, and 736–757 m for the rhyolite section and 955–976, 1178–1189, 1257–1273, and 1287–1356 m for the sandstone/mudstone section (Figure 6b,c). The highest temperature was 69.51 °C, at a depth of 1502 m. The thermal gradient of BH–2 was 22.89 °C/km at a depth of 130–370 m, 26.78 °C/km at a depth of 370–538 m, 64.52 °C/km at a depth of 538–770 m, and 33.33–35.67 °C/km at a depth of 770–1502 m.

The EC values at a depth of 122–1490 m ranged from 801 to 3758 μS/cm. The most distinctive anomalous zones of the EC were found at depths of 578–581, 622–627, 665, 737, 762, 943–955, 966–972, 1170, 1181, 1251, 1273, 1284, and 1360 m (Figure 6d).

The semiconsolidated mudstone had an average natural gamma radioactivity of 51 American Petroleum Institute (API) units at a depth of 0–356.7 m; lapilli tuff had 71 API units at a depth of 356.7–423.7 m; tuffaceous sediment had 103 API units at a depth of 423.7–426.4 m and 70 API units at a depth of 450.6–461.8 m; crystal tuff had 93 API units at a depth of 426.4–450.6 m; sandstone/mudstone had 72 APIs unit at a depth of 461.8–521.6 m and 56 API units at a depth of 876.4–1501.3 m; rhyolite had 64 API units at a depth of 521.6–701.1 m, 63 API units at a depth of 708.0–744.5 m, 65 API units at a depth of 775.7–827.9 m, and 66 API units at a depth of 830.6–876.4 m; and, finally, basic dyke had 33 API units at a depth of 701.1–708.0 m, 35 API units at a depth of 744.5–775.7 m, and 36 API units at a depth of 827.9–830.6 m (Figure 6e). The lower natural gamma-ray intensities represented sandstone at a depth of 1200–1440 m based on the core samples, as well as the basic dyke at a depth of 1440–1490 m. The zone of 1200–1440 m approximately corresponded with an anomaly depth of 1245 m based on the caliper log. The trend of the caliper log was different from that for the EC log, which could not distinguish between sandstone and basic dyke at depths greater than 1200 m. Rhyolite at a depth of 520–825 m exhibited a distinct natural gamma-ray intensity. Three basic dykes intruding into the rhyolite could be easily distinguished from the rhyolite, by their low natural gamma-ray values. At a depth of 460–520 m, the natural gamma-ray values of the semiconsolidated mudstone overlying the rhyolite were approximately the same with those of sandstone.

4.2.3. BH-3

Anomalous zones of the caliper log appeared at depths of 399–408 and 455–468 m for the basic dyke; 531–545, 564–578, and 605–607 m for sandstone/mudstone; and 862–877 m for rhyolite (Figure 7a).

The temperature log zones were identified at a depth of 138–381 m in the semiconsolidated mudstone, tuffs (lapilli tuff, crystal tuff, and tuffaceous sediment), and basic dyke sections (Figure 7b,c). The highest temperature of BH-3 was 44.33 °C at a depth of 912.27 m. The thermal gradient of BH-3 was 51.51 °C/km at a depth of 138–381 m and 40.85 °C/km at a depth of 381–900 m. The temperature log data for BH-3 had low accuracy due to the effect of the drilling mud, as they were obtained during drilling.

The EC values at a depth of 118–911 m ranged from 969 to 1561 μS/cm. The most distinctive anomalous zones of the EC can be observed at depths of 134–140 and 381–400 m (Figure 7d).

The semiconsolidated mudstone had an average natural gamma radioactivity of 64 API units at a depth of 0–206.1 m; crystal tuff had 104 API units at a depth of 206.1–257.0 m, 127 API units at a depth of 292.6–376.1 m, and 133 API units at a depth of 409.3–428.3 m; tuffaceous sediment had 119 API units at a depth of 257.0–273.2 m; basic dyke had 33 API units at a depth of 273.2–292.6 m and 30 API units at a depth of 428.3–474.7 m; lapilli tuff had 100 API units at a depth of 474.7–524.1 m; sandstone/mudstone had 77 API units at a depth of 524.1–628.7 m and 140 API units at a depth of 901.1–911.2 m; and, finally, rhyolite had 90 API units at a depth of 632.1–901.1 m (Figure 7e). According to the natural gamma-ray log data, BH-2 and BH-3 had similar lithologies.

4.2.4. BH-4

Anomalies of the caliper log occurred in the lower zones at depths of 1548, 1562, 1588–1599, 1733–1738, 1967, and 1988.7 m for the andesitic breccia (Figure 8a). At 828–1000 m, the caliper measurement was not executed, owing to the conditions of the borehole.

Abrupt zones of temperature changes were identified at depths of 303–350 m in the crystal tuff section, 394–480 m in the lapilli tuff section, and 567–600 and 873–967 m in the sandstone/mudstone section (Figure 8b,c). The thermal gradient of BH-4 was 20.92 °C/km at a depth of 115–225 m, 48.78 °C/km at a depth of 225–490 m, and 35.71 °C/km at a depth of 490–1475 m.

The EC values ranged from 513 to 2999 μS/cm between 106 and 1482 m. The most distinctive anomalous zones of the EC were found at depths of 301, 367, 394–397, 416, 585–598, 850, 967–970, 1383, and 1468 m (Figure 8d).

The semiconsolidated mudstone had an average natural gamma radioactivity of 24 API units at a depth of 0–217.4 m; crystal tuff had 133 API units at a depth of 217.4–252.4 m and 137 API units at a depth of 297.3–359.5 m; sandstone/mudstone had 37 API units at a depth of 443.3–480.1 m and 78 API units at a depth of 547.7–1361 m; lapilli tuff had 102 API units at a depth of 270.3–279.3 m, 120 API units at a depth of 368.6–443.3 m, and 95 API units at a depth of 480.1–547.7 m; and, finally, andesitic breccia had 48 API units at a depth of 1361–1.996 m (Figure 8e). The low natural gamma-ray intensity at a depth of 740–810 m indicated the occurrence of sandstone, whereas the higher values at 1285 m indicated shale and andesitic breccia, which do not exist in the other boreholes; a high intensity was observed at a depth of 1360 m for BH-4.

4.3. Groundwater Flow Characteristics

The groundwater flow based on the temperature log data was interpreted using the 1D steady-state heat transfer solution [11,27] and by matching the observed temperature data with appropriate types of curves, considering the lithology log and anomalies of the caliper and EC logs in the process. By matching to the most appropriate type of curve, the corresponding β value, groundwater flow velocity, and upward/downward flow directions were determined with the residual sum of squares (RSS), 0.008–2.013 (Table 1). The thermal conductivity values of 2.16–3.95 W/mK (Table 1) and thermal diffusivity values of 0.71–1.60 mm²/s were determined from 42 samples. One sample was taken from the basic dyke, seven samples were taken from rhyolite, and thirty-four samples were taken from sandstone/mudstone layers. The unmeasured thermal conductivity values of the dacite and crystal tuff were based on representative values in South Korea (Table 1) [52]. The determined thermal diffusivity and thermal conductivity values seem reasonable by the comparison between the reference values [53,54,55] (

Table 2. Comparison of reference and measured values of thermal diffusivity and thermal conductivity for different rock types.

Rock Types	Thermal Diffusivity (mm2/s)			Thermal Conductivity (W/mK)
	[53,54]	[55]	Measured Data	[53,54]	[55]	Measured Data
Dacite	-	1.06	-	-	2.73	-
Sandstone/ mudstone	0.80–0.96	0.80–1.51	1.22	2.04–2.47	2.04–3.86	3.01
Basic dyke	-	-	0.74	-	-	2.86
Rhyolite	1.40	1.20	1.51	3.59	3.08	3.87
Tuff	0.80	1.17	-	2.05	3.00	-

).

Williams [56] reported that a harmonic mean of measured thermal conductivities could appropriately consider the anisotropy of stratified formations. In this study, the harmonic means of the thermal conductivity values were applied for three or more samples of the rocks. A specific heat of water in the range of 4.18–4.19 kJ/kgK and a density in the range of 980–998 kg/m³ at 22–65 °C were applied for the 1D steady-state heat transfer.

In BH-1, upward flows occurred along the fractures with dip angles of −60° to 70° in four sections: 358–410 m along a fracture zone of 410.3–410.7 m in the dacite and tuffaceous sediment section with a flow velocity of −2.27 × 10⁻⁵ m/s; 488–514 m along a fracture zone of 511–514.8 m in the sandstone/mudstone with a flow velocity of −5.93 × 10⁻⁵ m/s; 540–566 m along a fracture zone of 523.8–569.5 m in the rhyolite, with a flow velocity of −7.51 × 10⁻⁵ m/s; and 890–998 m along a fracture zone of 977.2–979.2 m in the sandstone/mudstone, with a flow velocity of −7.86 × 10⁻⁶ m/s. Alternatively, a downward flow occurred at 586–646 m along a fault zone of 572.6–572.9 m in the rhyolite, with a dip angle of −50° and flow velocity of 1.09 × 10⁻⁴ m/s.

In BH-2, upward flows with dip angles of −10°–60° occurred in four sections: 618–652 m along a highly fractured zone of 525–900 m in the rhyolite with a flow velocity of −5.16 × 10⁻⁴ m/s, corresponding to the EC anomaly at 655 m; 736–757 m along a highly fractured zone of 525–900 m in the basic dyke with a flow velocity of −1.41 × 10⁻⁴ m/s, corresponding to the EC anomaly at 762 m; 1178–1189 m along a highly fractured zone of 1155–1165 m in the sandstone/mudstone with a flow velocity of −5.90 × 10⁻⁴ m/s, corresponding to the EC anomaly at 1181 m; and 1287–1356 m along the fault of 1310 m in the sandstone/mudstone with a flow velocity of −1.97 × 10⁻⁵ m/s, corresponding to the EC anomaly at 1360 m. Downward flows with dip angles of −17° to 28° occurred in three sections: 580–610 m along the highly fractured zone of 525–900 m in the rhyolite with a flow velocity of 3.16 × 10⁻⁴ m/s, corresponding to the EC anomaly at 578–581 m; 955–976 m along the bedding plane in the sandstone/mudstone with a flow velocity of 2.39 × 10⁻⁴ m/s, corresponding to the EC anomaly at 943–955 m in the sandstone/mudstone; and 1257–1273 m along a fault of 1240 m in the sandstone/mudstone with a flow velocity of 4.69 × 10⁻⁴ m/s, corresponding to the EC anomaly at 1251 m. Figure 9 shows downward flow in the zone of 1257–1273 m with β = 3.90 (RSS = 0.031). The deviations between the type of curve and measured temperature values, especially in the section of 0.2 < z_D < 0.4 in Figure 9 can be explained by the deviations of the centering of the loggers [57] or change of fracture direction or interconnection with other facture sets [27]. The deviation can also be noticed on the EC profile in Figure 6.

In BH-3, an upward flow with a dip angle of −90° occurred at 138–381 m along a fault of 400–403 m in the basic dyke with a flow velocity of −5.27 × 10⁻⁶ m/s, corresponding to the EC anomaly at 381–400 m.

In BH-4, an upward flow with a dip angle of approximately 65° occurred 567–600 m along a fault of 501–502 m in the sandstone/mudstone with a flow velocity of −8.07 × 10⁻⁵ m/s, corresponding to the EC anomaly at 585–598 m. Alternatively, downward flows with dip angles of 30°–90° occurred in three sections: 303–350 m along the fault of 281.9–282 m in the crystal tuff with a flow velocity of 2.17 × 10⁻⁴ m/s, corresponding to the EC anomaly at 301 m; 394–480 m along the fault of 402–402.1 m in the crystal tuff with a flow velocity of 7.96 × 10⁻⁵ m/s, corresponding to the EC anomaly at 394–397 m; and 873–967 m along the fracture zone of 1000–1127 m in the sandstone/mudstone with a flow velocity of 4.28 × 10⁻⁵ m/s, corresponding to the EC anomaly at 850 m. Figure 9 shows the upward flow in the zone of 567–600 m, with β = −3.39 (RSS = 0.010).

5. Discussion

Conventional geophysical logs can qualitatively indicate the fractures and alteration zones intersecting the boreholes [58,59]. In this study, we characterized deep groundwater flow and thermal properties by using geophysical logs and lithological logs of boreholes. The fracture zones and thermal properties were identified based on the geophysical logs, following the research of Paillet [60,61] and Katsube and Hume [62].

The geophysical logs and lithology profiles indicated that groundwater flows into, and out of, the boreholes in the basic dyke, rhyolite, and sandstone/mudstone that contain highly fractured zones. The caliper log identified the anomaly zones associated with highly fractured zones in the basic dyke, rhyolite, sandstone/mudstone, and andesitic breccia. The results of the strata analysis using natural gamma-ray log and lithology log suggest that the geological strata from the upper Mesozoic sedimentary rocks to the upper part included semiconsolidated mudstone, tuffs (lapilli tuff, crystal tuff, and tuffaceous sediment), rhyolite, and sandstone/mudstone; moreover, the strata possessed good horizontal connectivity, with possible fracture connections in the rock around the boreholes. The depth of the Cretaceous semiconsolidated mudstone layer increased in depth from north (BH-3 and BH-4) to south (BH-1 and BH-2). The distribution of the tuffaceous sediment decreased and that of the crystal tuff increased from south to north.

By comparing the caliper log and lithology data, anomalies appeared for BH-2 at the following depths: 424, 461, and 504 m for the rhyolite; 704, 731, and 760 m for the basic dyke; and 1087, 1125, and 1245 m for the sandstone/mudstone. This is similar to the anomalous zone of 1200–1440 m obtained based on the natural gamma-ray log. The anomalies in BH-3 appeared at a depth of 403 m for the basic dyke, 620 m for the sandstone/mudstone, and 862 and 877 m for the rhyolite. In contrast, the anomalies in BH-4 occurred at greater depths, specifically, at 1589, 1738, and 1989 m for the andesitic breccia. In BH-2, the three basic dykes (intruding into the rhyolite) were distinguished easily from the rhyolite, owing to their low natural gamma-ray intensities.

By means of the lithological log and geophysical logs, several aquifer zones were identified, at depths of 350–410, 510–590, 600–650, 730–760, 870–890, 950–1000, and 1150–1360 m, as well as hydraulically active fractures intersecting the deep boreholes. An anomaly in the EC values can represent an inflow or outflow point of fresh and/or saline water in the borehole [45].

Uncertainty can be induced in the upward or downward flows of the 1D steady-state heat transfer model from the deviations in the centering of the loggers [57] or change of fracture direction or interconnection with other facture sets [27] in the boreholes, e.g., BH-2 in Figure 9. The average angle of the fracture zone, $\omega$, representing the average direction of groundwater flow along the fracture, is not known a priori [27]. Uncertainty also exists in the groundwater flow owing to indistinct nonlinear temperature zones; these can be more clearly identified through a detailed investigation of the fractures and faults in deep boreholes.

6. Conclusions

The lithological and geophysical logs of the four deep boreholes could identify the deep groundwater flows in the Pohang area that experienced the youngest tectonic movements in the Korean peninsula. The groundwater flow directions and β values of the boreholes were estimated by using a 1D steady-state heat transfer model. The geological stratigraphy of the four boreholes primarily contained Tertiary semiconsolidated mudstone, dacite, tuffs (crystal tuff, lapilli tuff, and tuffaceous sediment), basic dyke, sandstone/mudstone, andesitic breccia, and granodiorite. BH-1, BH-2, and BH-3 had similar stratigraphy overall, and no rhyolite was found in BH-4 at depths greater than 600 m. The semiconsolidated mudstone was located at a greater depth in BH-1 and BH-2 than in BH-3 and BH-4.

In the interpretation of the groundwater flow, the fractures and faults induced nonlinear temperature changes and upward and downward flows. Moreover, several aquifer zones were identified at depths of 350–410, 510–590, 600–650, 730–760, 870–890, 950–1000, and 1150–1360 m, as well as hydraulically active fractures intersecting the deep boreholes. In BH-1, upward flows occurred at four depths, at 358–410, 488–514, 540–566, and 890–998 m, whereas a downward flow occurred at 586–646 m. Similarly, in BH-2, upward flows occurred at four depths, at 618–652, 736–757, 1178–1189, and 1287–1356 m, whereas a downward flow occurred at 580–610 and 955–976 m. In BH-3, an upward flow occurred at 138–381 m, with a linear temperature change. In BH-4, an upward flow occurred at 567–600 m, whereas downward flows occurred at 303–350, 394–480, and 873–967 m.

In conclusion, the upward and downward flows were identified by solving the 1D steady-state heat transfer equation, in combination with the lithological log and geophysical logs. However, owing to the loss of the core recovery, slight discrepancies were observed between the geophysical logs and lithology log. The groundwater flow occurred mostly along the fault and fracture zones in the rhyolite, basic dyke, and sandstone/mudstone. In addition, by using the lithological log and geophysical logs, several aquifer zones and hydraulically active fractures intersecting the deep boreholes were identified.

Upward or downward flows based on the 1D steady-state heat transfer model can be imprecisely identified, owing to the deviations in the centering of the loggers in the boreholes. There also exists uncertainty regarding groundwater flow along depth, since more complex groundwater flow patterns can be detected through detailed investigations of fractures and faults, and these could not be considered in the deep boreholes. Despite the uncertainty, analyses combining geophysical logs and a lithological log, as well as a 1D steady-state heat transfer model, can contribute in understanding the deep fractures/fault distributions/thermal and flow characteristics of boreholes in relation to hydrogeological activities, through the analysis of the deep fractures/faults and geothermal potential around the study area.