Hydrochemical and Isotopic Difference of Spring Water Depending on Flow Type in a Stratigraphically Complex Karst Area of South Korea

Characterizing the subsurface flow in karstic areas is challenging due to distinct flow paths coexisting, and lithologic heterogeneity makes it more difficult. A combined use of hydrochemical, environmental isotopic, and hydrograph separation study was performed to understand the subsurface flow in a karst terrain where Ordovician carbonate rocks overlie Jurassic sandstone and shale along thrusts. Spring water collected was divided into Type Ⅰ (n = 11) and Ⅱ (n = 30) based on flow patterns (i.e., low and high discharge, respectively). In addition, groundwater (n = 20) was examined for comparison. Three Type Ⅱ springs were additionally collected during a storm event to construct hydrographs using δ¹⁸O and δD. As a result, Type Ⅱ had higher electrical conductivity, Mg²⁺, HCO₃⁻, and Ca²⁺/(Na⁺ + K⁺) than Type Ⅰ and was mostly saturated with calcite, similar to deep groundwater. The hydrochemical difference between Types Ⅰ and Ⅱ was opposite to the expectation that Type Ⅱ would be undersaturated given fast flow and small storage, which could be explained by the distinct geology and water sources. Most Type Ⅱ springs and deep groundwater occurred in carbonate rocks, whereas Type Ⅰ and shallow groundwater occurred in various geological settings. The carbonate rocks seemed to provide conduit flow paths for Type Ⅱ given high solubility and faults, resulting in 1) relatively high tritium and NO₃⁻ and Cl⁻via short-circuiting flow paths and 2) the similar hydrochemistry and δ¹⁸O and δD to deep groundwater via upwelling from deep flow paths. The deep groundwater contributed to 83–87% of the discharge at three Type Ⅱ springs in the dry season. In contrast, Type Ⅰ showed low Ca²⁺ + Mg²⁺ and Ca²⁺/(Na⁺ + K⁺) discharging diffuse sources passing through shallow depths in a matrix in mountain areas. Delayed responses to rainfall and the increased concentrations of contaminants (e.g., NO₃⁻) during a typhoon at Type Ⅱ implied storage in the vadose zone. This study shows that hydrochemical and isotopic investigations are effective to characterize flow paths, when combined with hydrograph separation because the heterogenous geology affects both flow paths and the hydrochemistry of spring water passing through each pathway.

Introduction

Karst aquifers differ from unconsolidated granular aquifers regarding subsurface flow behavior because of various types of porosity, e.g., matrix, fractures, and conduits (Ghasemizadeh et al., 2012; Hartmann et al., 2014; Demiroglu, 2016). In addition, the epikarst plays a significant role in karst recharge and storage (Aquilina et al., 2005; Aquilina et al., 2006; Doctor et al., 2006; Chang et al., 2019). Moreover, lithologic heterogeneity often makes the flow system in karst aquifers complicated (Mocior et al., 2015; Lorette et al., 2018; Frank et al., 2019; Zhu et al., 2020).

In the karst area, groundwater recharge occurs by either diffuse infiltration through the fissured matrix, called diffuse-dominated, or direct point infiltration into the conduit network, called conduit-dominant. The diffuse-dominated flow is associated with slow and delayed responses of water seepage from the aquifer matrix, less connected fractures, or overlying soils, whereas the conduit-dominated flow is fast through the interconnected voids and passages of the aquifer (Schilling and Helmers, 2008). Theoretically, the diffuse laminar flow is distinguished from the turbulent conduit flow based on flow velocities (e.g., a few centimeters or meters per day for diffuse flow in Atkinson, 1977; a few m/h to several hundred m/h for conduits in Rusjan et al., 2019) or Reynolds numbers (e.g., less than 2,000 for laminar flow in Grasso et al., 2003; Ghasemizadeh et al., 2012; over 100 for conduit flow in Lee et al., 2006). The diffuse flow in general is assumed to undergo water-rock interactions, increasing the contents of dissolved ions (e.g., Adji et al., 2016), whereas the conduit-dominant flow is habitually undersaturated with respect to calcite (e.g., De Rooij and Graham, 2017). Of them, conduit flow makes the karst aquifer vulnerable to contamination through direct and fast infiltration (Stueber and Criss, 2005; Ghasemizadeh et al., 2012; Parise et al., 2015) or to rock collapse due to open caves (Waltham and Lu, 2007; Pogačnik et al., 2017). Thus, an understanding of conduit flow paths is essential when managing water resources or building facilities (e.g., dams) in the karst area (Milanović, 2021).

Subsurface flow in karst aquifers is mainly characterized by monitoring the hydrological (e.g., discharge), hydrochemical (e.g., cations and anions), and isotopic (mainly, δD, δ¹⁸O, and δ¹³C) variability of spring water that is groundwater naturally emerging to the surface (e.g., Lee and Krothe, 2001; Aquilina et al., 2005; Barbieri et al., 2005; Aquilina et al., 2006; Doctor et al., 2006; Moore et al., 2009; Jiang et al., 2013; Lorette et al., 2018; Gil-Márquez et al., 2019; Rusjan et al., 2019; Palcsu et al., 2021). For instance, Moore et al. (2009) identified water sources (e.g., deep-water upwelling) and their mixing and interactions with aquifer materials based on hydrochemical data to understand the relationship between spring characteristics and the complexities of karst aquifers. Doctor et al. (2006) used the δD and δ¹⁸O of water and δ¹³C of dissolved inorganic carbon (DIC) as well as major ion chemistry to estimate mixing proportions among multiple sources, and found that water released from storage within the epikarst may comprise as much as two-thirds of conduit flow in a karst aquifer following rainfall. Palcsu et al. (2021) found the 10% of the discharge water with a residence time of half a year and the contribution of an older component based on the seasonality of δ¹⁸O in 2015–2018 at a karst spring with an average discharge of 20 L/min in Hungary. Lorette et al. (2018) found high Mg²⁺ and low dissolved oxygen (DO) during low-water periods and explained the uncommon hydrochemistry using the mixture with a Jurassic confined aquifer in western France.

In those studies, geochemical modelling and saturation index (SI) computation were often applied to interpret the monitoring data (e.g., Barbieri et al., 2005; Peyraube et al., 2012) because the source of acidity for carbonate dissolution and the geochemical conditions for geochemical facies as well as the impact of human activities can be evaluated by examining the relationships between solution concentrations, partial pressure of CO₂ (P_CO2), and SI of calcite (SI_calcite). Also, principal component analysis (PCA) was extensively used to infer hydrogeochemical processes or water sources (e.g., Doctor et al., 2006; Moore et al., 2009; Mudarra et al., 2012; Lorette et al., 2018; Gil-Márquez et al., 2019). Mudarra et al. (2012) included P_CO2 and SI_calcite in PCA.

Another method to understand the groundwater flow in a karst area is an examination of spring responses to rainfall to quantify the amount of rainwater (or baseflow) contributing to the discharge. Hydrograph separation is performed using various data including hydrochemical parameters (e.g., NO₃⁻, Cl⁻, and specific conductance by Schilling and Helmers (2008)) and stable isotopes (e.g., δD and δ¹⁸O by Aquilina et al. (2005), Klaus and MacDonnell (2013)) to apportion the discharge among multiple water sources. Vesper and White (2004) noted that karst spring hydrographs reflect not only the changing mix of baseflow and storm flow, but also a shift in recharge sources during a storm event. Grasso et al. (2003) coupled the analysis of karst spring hydrographs and chemographs and evaluated the geometric dimensions of submerged karstic networks.

We noted that the combined use of hydrochemical and isotopic compositions and hydrograph separation with the help of PCA would be effective to describe the subsurface flow in karstic areas as in Gil-Márquez et al. (2019), who used hydrodynamics, hydrochemistry, and environmental isotopes to devise a hydrogeological conceptual model of evaporite-karst springs. In particular, the multi-approach is expected to help delineate flow paths when the karstic flow system is intricate with lithologic heterogeneity. Besides, spatial studies were rarely conducted for the hydrochemical and isotopic differences of spring water depending on flow types (e.g., low and high discharge), whereas the temporal variation of discharge and hydrochemistry as a function of the hydrological regime was performed by many researchers, including Bicalho et al. (2012)Adji et al. (2016), Lorette et al. (2018), Gil-Márquez et al. (2019), Frank et al. (2019), and de Filippi et al. (2021).

Therefore, we conducted hydrochemical, environmental isotopic, and hydrograph separation studies to assess the subsurface flow in a stratigraphically complex karst area with thrusts where both low discharge and high discharge springs occurred (Figure 1). It was assumed that the low discharge and high discharge come through diffuse and conduit pathways, respectively, because the flow in karst aquifers is turbulent in conduits such as caves and sinkholes, while diffusive through the fissured matrix or intergranular unfractured bedrock (Demiroglu, 2016; Chang et al., 2019). Specifically, springs were divided into Type Ⅰ and Ⅱ when sampled based on the manner by which the discharge occurred (i.e., low and high, respectively) as in Mocior et al. (2015) who surveyed the spatial distribution of springs in mountains and categorized them by discharge. Then the hydrochemical and isotopic properties of spring water were characterized depending on the flow type and compared with those of groundwater to understand the flow path. Then baseflow separation was conducted for three Type Ⅱ springs by monitoring δ¹⁸O and δD data during a storm event to estimate the proportion of baseflow in the conduit flow-dominated spring discharge. Lastly, the subsurface flow and the storage of water and solutes in the vadose zone were discussed.

FIGURE 1

FIGURE 1. Geologic map of the study area in South Korea: (A) plan view; (B) cross sections along A-A’, B-B’, C-C’ and D-D’ in (A). In (A), the white area indicates Kosong limestone (Ogl) of unknown age. Sampling locations are also shown, while the names are labeled only for the samples mentioned in the text. All sample names are shown on the land use map in Supplementary Figure S2. Springs were divided into Type Ⅰ and Ⅱ when sampled based on the manner by which the discharge occurred (i.e., low and high, respectively).

The contribution of this work is the characterization of groundwater flow systems in a stratigraphically complex karst aquifer with thrusts using the spatial distribution of hydrochemical and isotopic data with statistical and hydrographic separation analysis in South Korea which has few studies on karstic aquifers. It should be noted that a multipurpose dam was planned to be built on the Dong River on the Bansong Formation (Jbs and Jbc; Figure 1) in the early 1990s for water supply, flood management, and hydropower generation, but the plan was canceled in 2000 because of environmental issues. The field work for this study was conducted during the feasibility study of the dam since karst structures in heterogeneous lithological units may cause geohazards for civil infrastructures (Pogačnik et al., 2017). The data has not been properly interpreted until this study because of social sensitivity. Although the data are >20 years, the data interpretation is expected to provide information on the subsurface flow in this stratigraphically complex karst aquifer given the geologic time scale and no significant land use change over 20 years in the area (Supplementary Figure S1).

Study Area

Geography and Climate

The study area is located in the Yeongwol-Danyang basin at the middle eastern part of South Korea (Figure 1). The topographic elevation ranges from 200 to 800 m above sea level (asl). In terms of geomorphology, high and steep mountains are located within the study area, and flat terranes are developed on a mountain ridge at high altitude (Supplementary Figure S2), indicating repeated cycles of uplifting and erosion in this site over geological time. The Dong River meanders along narrow and deep valleys formed between high and steep mountains (Figure 1A). The land along the river has been used for agriculture since the late 1980s (Supplementary Figure S1), which implies that the study area was exposed to agrochemicals when water samples were collected as it is today.

The study area has four distinct seasons. The summer in June to August is hot and humid, whereas the winter in December to February is cold and dry. In summer, typhoons bring strong winds and heavy rain. Specifically, in 1999, the highest monthly average temperature was 23.6°C in July and August, while the lowest was −2.7°C in January. Similarly, the highest and lowest in 2000 were 24.4°C in July and −3.3°C in February, respectively. The total annual precipitation was 1,360.6 mm in 1999 and 1,060.8 mm in 2000, of which 70 and 78% occurred in June to September in 1999 and 2000, respectively. The precipitation during the study period is shown in Figure 2. The total amount of rainfall was 115 mm over approximately 24 h and the maximum rainfall intensity was <10 mm/h during a typhoon called Saomai in 2000 (Figure 2B).

FIGURE 2

FIGURE 2. Hydrological responses of three Type Ⅱ springs from March to July in 1999 (A) and during the typhoon event named Saomai in 2000 (B). The blue bar indicates precipitation.

Geology

The geology of the study area consists of the Cambro-Ordovician Joseon Supergroup (mainly the Maggol (Omg) and Yeongheung (Oy) Formations), Carboniferous to Permian Pyeongan Supergroup (e.g., Hongjeom series (Ch)), and Late Triassic to Early Jurassic Daedong Supergroup (Bansong Formation (Jbs and Jbc)) as in Figure 1 (Han et al., 2006). The age of Kosong limestone (Ogl) on the east (the white area in Figure 1A) is unknown, and had been suggested to be Ordovician or later (Seo, 1997).

The Maggol Formation (Omg) comprises a thick (250–400 m) sequence of carbonate rocks with a variety of lithology such as calcareous mudstone, dolostone, lime conglomerate, lime breccia, bioclastic grainstone, and oolitic grainstone, while the Yeongheung Formation (Oy; ca. 400 m thick) consists of massive to thick-bedded, light to dark gray dolostone in the lower part and bluish gray limestone in the upper part (Chough et al., 2000). These carbonate rocks occupy most of the study area (Figure 1A), and are characterized by fractures (Figure 1B) and sinkholes, caverns, and dolines (Supplementary Figure S2). The carbonate rock aquifers yield a usable quantity of water in the study area (http://www.gims.go.kr). Park et al. (2011) reported the hydraulic conductivity of 0.004–1.1 m/day (mean of 0.078 m/day) in carbonate aquifers near the study area. Recently a tracer test was performed using uranine at the west of the study area (out of the map in Figure 1A) in January 2021 by a private company which plans to build a waste landfill in a historical limestone mine site. The uranine leaked to a small stream nearby 3 days after being released at the abandoned mine site (https://news.kbs.co.kr/news/view.do?ncd=5158987), implying hydraulic connectivity.

Meanwhile, the Bansong Formation (Jbs and Jbc) consists of a conglomerate succession at the base and metamorphosed shale and sandstone with thin coal seams at the upper part, and unconformably overlies the carbonate sequence of the Joseon Supergroup or the siliciclastic sequence of the Pyeongan Supergroup (Han et al., 2006). The Bansong Group is highly metamorphosed and deformed in places such as the northeast-trending regional folds and thrusts developed during the Daebo tectonic event (Han et al., 2006; Ree et al., 2009), and represents syntectonic sedimentation during the late Early Jurassic in a compressional regime. Specifically, the sequence is structurally overlain by the Maggol Formation (Omg) along the Gongsuweon thrust and by the Yeongheung Formation (Oy) along the Deokpori thrust (Figure 1B). According to Ree et al. (2009), the deposition of the Bansong Group and the thrusting were coeval during the middle Mesozoic Daebo tectonic event.

Given that springs in general originate along geologic contacts, faults, and ground depression (Fiorillo et al., 2018) and the layering of rocks relative to the slope of mountains affect both the number and discharge of springs in mountains (Mocior et al., 2015), the overturned strata (e.g., Cambro-Ordovician Yeongheung Formation (Oy) overlying Jurassic Bansong Formation (Jbs) along the thrust in Figure 1B) in this mountainous area can provide a pathway for deep groundwater upwelling as in Lorette et al. (2018) who showed deep water rising at springs through the faulted anticline structure of Périgueux, France. Similarly, Zhu et al. (2020) observed ascending springs along secondary tectonic fissures in the contact zone of limestone and magmatic rocks invaded as well as contact springs where the permeability of the local contact zone was enhanced. According to Zhu et al. (2020), the two contact springs with relatively large and stable flow rates reflected the high proportion of karst water replenishment from the deep circulation.

Methods

Sampling

Water samples were collected at two steps. First, various kinds of water were collected for major chemical compositions in the dry season (May/June; n = 44) and wet season (July; n = 27) in 1999 (Figure 1), including rainwater (n = 2 from RW1 in Figure 1A by installing a rainfall collector), spring water (n = 41 from 30 sites), surface water (n = 8 from two sites at Dong River and four sites at streams), shallow groundwater (n = 11 from seven preexisting wells), and deep groundwater (n = 9 from six preexisting wells) to assess the hydrochemistry of each water group. Sampling dates were summarized in Supplementary Table S1 where sampling sites visited in both dry and wet seasons were in bold. In addition, water samples were collected for δD and δ¹⁸O of water (n = 18) in the dry season and for δ¹³C of DIC (δ¹³C_DIC; n = 10) and tritium (n = 11) in the wet season to identify water sources. The lithology of each sampling site was assessed based on a geologic map (Figure 1), although the surface geology might not represent the subsurface geology through which the water passed, in particular in this stratigraphically complex study area with thrusts.

Spring water samples were divided into Type Ⅰ (n = 11; much lower than the rates in Figure 2A) and Type Ⅱ (n = 30; discharge rates similar to those in Figure 2A (> 480 L/min)) based on the flow type observed during sampling to assess the hydrochemical and isotopic differences of spring water depending on flow patterns and to determine the factors contributing to the differences in the study area with stratigraphic overturning along thrusts. Type Ⅱ springs were expected to exhibit large variations in discharge and chemical composition through time (Moore et al., 2009). In fact, Type Ⅱ quickly responded to rain events as in Figure 2B, and the discharge rate was over 1.7 m³/s at S46 (Type Ⅱ) during the storm event in 2000. However, Reynolds numbers were not calculated to confirm the laminar or turbulent flow due to little information on conduit sizes, and discharge rates were measured only at three Type Ⅱ springs in Figure 2 due to the limited accessibility to the springs in the wild areas, which is the limitation of this study.

In the second stage, a field investigation was conducted from September 15 to 19 in 2000 to characterize the short-term variability of isotopic compositions over the course of the typhoon event named Saomai. Three Type Ⅱ springs were chosen (S41, S45, and S46 in Figure 1A) from the Yeongheung (Oy) Formation to assess water sources through conduits. They were selected because of relatively easy accessibility and discharge measurement. Spring water samples for stable isotopes (δD and δ¹⁸O) were collected two or three times a day from the three springs (a total of 10 samples for each spring), and discharge rates were obtained through a constant cross section using a flowmeter (General Oceanics). In addition, a rainwater sample was collected at a relatively high elevation point (RW1 in Figure 1A).

Water samples were collected in 60 ml Nalgene^® bottles for the analysis of major compositions (Na⁺, K⁺, Ca²⁺, Mg²⁺, SO₄²⁻, NO₃⁻, Cl⁻, F⁻, and SiO₂) and δD and δ¹⁸O, while 1 L of water was collected for tritium analysis. Water samples were filtered through 0.45 μm pore size cellulose nitrate membrane filters (Whatman^®) to eliminate suspended materials and then stored in the polyethylene bottles. To the samples for cation analysis, concentrated nitric acid was added to keep pH <2. All water samples were stored at a constant low temperature of 4°C until analysis in the laboratory. In addition, DIC was obtained for δ¹³C_DIC analysis by precipitating barium carbonate (BaCO₃) with adding sodium hydroxide (NaOH) powder and barium chloride (BaCl₂) liquid to a 1 L water sample.

Analysis

In-situ measurements, including temperature, pH, redox potential (Eh), DO, and electrical conductivity (EC), were determined using a pH/Eh meter (Model Orion 290A), DO meter (Model Orion 835), and conductivity meter (Model Orion 130) in the field, while alkalinity was determined by titration in the field and re-examined with the Gran method in the laboratory.

Major cations and SiO₂ and anions were analyzed by inductively coupled plasma atomic emission spectrometry (ICP-AES; Model Perkin Elmer Optima 3000XL) and ion chromatography (IC; Model Dionex 120), respectively, at the Center for Mineral Resource Research (CMR) in Korea University. The δ¹⁸O and δD of water and δ¹³C_DIC were determined with an isotope ratio mass spectrometer (Finnigan MAT 252) at CMR. The δ¹⁸O and δD values were reported related to the Vienna Standard Mean Ocean Water (V-SMOW) with a precision of ± 0.1 and ± 1‰, respectively, while δ¹³C_DIC was expressed related to the Pee Dee Belemnite (PDB) and its standard deviation was ± 0.05‰. The tritium (³H) contents in water samples were measured at the Korea Atomic Energy Research Institute with a liquid scintillation counter (Model Packard 3255) after a ca. 600 g water sample was electrolytically condensed to be 20 g. The analytical error was within approximately 0.2 TU.

Thermodynamic, Statistical, and Hydrographic Separation Analysis

Partial pressure of CO₂ (P_CO2) and the saturation index of calcite (SI_calcite) and dolomite (SI_dolomite) in water samples were calculated using the physicochemical data of water samples and the PHREEQC program (Parkhurst and Appelo, 2013), by assuming an equilibrium state. Then, the binary relations between SI_calcite and SI_dolomite (Langmuir, 1971), between SI_calcite and P_CO2 (Doctor et al., 2006; Lacelle et al., 2008; Frondini et al., 2019), and between pH and HCO₃⁻ (Langmuir, 1971; Lacelle et al., 2008) were drawn to assess the hydrochemical evolution with carbonate minerals. Groundwater flows at relatively high velocities in a karst area, and it is difficult to assume the equilibrium state. However, equilibrium thermodynamics provide quantitative information about geochemical conditions, including the open or closed system dissolution defined by Langmuir (1971). Once the equilibrium is interpreted, the hydrochemical and hydrodynamic evolution can be deduced.

Principal component analysis (PCA) was applied to the physicochemical parameters and SI_calcite and SI_dolomite of water samples (n = 69), excluding rainwater to assess major hydrogeochemical processes. PC loading and scores are considered together to assess major hydrochemical processes in the studied area, while water samples were plotted in the plot of PC scores to differentiate a major hydrochemical process for each water group (Moore et al., 2009; Bicalho et al., 2012; Lorette et al., 2018; Kim et al., 2021). PCA was conducted using the latest MATLAB software (https://www.mathworks.com) after log-transformation except pH.

Hydrograph separation curves were made for the three Type Ⅱ springs using δ¹⁸O and δD to better understand the discharge at Type Ⅱ springs. A two-component mixing model was applied with rainwater (R) and pre-storm water (PS), and then the contribution of R in a spring water sample (M) was evaluated as:

$Q_{\mathit{R}}=Q_M \times \left({\delta }_M-{\delta }_{\mathit{PS}}\right)/\left({\delta }_{\mathit{R}}-{\delta }_{\mathit{PS}}\right)$

where Q is discharge (L/s) and $\delta$ is the isotopic composition. Changes in $Q_M$ and $Q_{PS}$ were plotted with time to see the different proportion of $Q_{PS}$ in the total discharge during a storm event. Note that soil water and epikarst water were excluded due to little information despite their contributions to discharge (Lee and Krothe, 2001; Aquilina et al., 2005; Aquilina et al., 2006). Instead, the water and solutes stored in the vadose zone were discussed with the concentrations of anthropogenic contaminants (e.g., NO₃⁻) given the contaminant storage occurring in the vadose zone including soil and epikarst (Ghasemizadeh et al., 2012) as well as water storage (Perrin et al., 2003; Aquilina et al., 2006; Doctor et al., 2006; Jacob et al., 2008).

Results and Discussion

Hydrochemical and Isotopic Characteristics

Physicochemical data for each water group are summarized in Table 1. The water samples from both dry and wet seasons were combined to characterize the hydrochemistry because of no temporal variation between two seasons on the Piper’s diagram with NO₃⁻ on the Cl⁻ axis (Figure 3), which showed that most of the water samples had the Ca-(Mg)-HCO₃ type. Given that Types Ⅰ and Ⅱ showed different hydrochemical and isotopic compositions (Table 1 and Figures 3–6), the division seemed acceptable, although it was not confirmed using Reynolds numbers.

TABLE 1

TABLE 1. Summary of physicochemical data for various kinds of water groups from the Yeongweol area, South Korea.

FIGURE 3

FIGURE 3. Piper’s diagram with NO₃⁻ on the Cl⁻ axis showing the hydrochemical characteristics of various kinds of water groups collected during the dry season (A) and rainy season (B) in 1999. The letters on the diagram indicate the names of distinct samples.

Groundwater and Surface Water

The deep groundwater, whose average depth to water level and average water table were respectively 102.0 and 258.0 m asl, contained higher Mg²⁺ and HCO₃⁻ than the shallow groundwater with an average depth to water level of 11.3 m and an average water table of 265.6 m asl. Whereas the shallow groundwater contained higher Cl⁻ and NO₃⁻ than the deep groundwater (Table 1; Figure 4). The hydrochemical difference indicated substantial carbonate dissolution in the deep groundwater and the shallow groundwater vulnerable to contaminants from the ground surface because Mg²⁺ increases as the water-rock interaction proceeds in dolomitic geologic settings (e.g., Lorette et al., 2018; De Filippi et al., 2021), while Cl⁻ and NO₃⁻ can be indicators of agricultural contaminants (Schilling and Helmers, 2008; Kim et al., 2021). Consistently, the deep groundwater had lower tritium contents (3.3–6.7 TU) than the shallow groundwater (8.7–9.4 TU) (Figure 5A), indicating relatively old ages (Palcsu et al., 2021), and high δ¹³C_DIC up to −10.0‰ (Figure 5B), implying a larger carbonate dissolution effect. According to Marfia et al. (2004), the DIC in groundwater evolves initially under open system conditions (δ¹³C_DIC of −25 to −17‰), and later under closed system conditions as the δ¹³C_DIC increases (δ¹³C_DIC of –17 to –13‰) as a result of carbonate dissolution.

FIGURE 4

FIGURE 4. Comparison of concentrations: (A) Mg²⁺ vs. HCO₃⁻; (B) Ca²⁺ vs. EC; (C) NO₃⁻ vs. Cl⁻. In (A) and (C), the letters indicate the name of distinct samples.

FIGURE 5

FIGURE 5. Environmental isotropic compositions: (A) Tritium concentrations (n = 11); (B) carbon isotopes of dissolved inorganic carbon (δ¹³C_DIC; n = 10) collected during wet season. Hydrogen (δD) and oxygen isotopes (δ¹⁸O) of water collected in the dry season in 1999 (n = 18) (C) and of three Type Ⅱ springs collected during the storm event in 2000 (D). In (D), the Groundwater (average) came from deep groundwater in (C) (δD = −70.2‰; δ¹⁸O = −9.6‰; n = 3) and was used to represent the pre-storm water (PS) for hydrograph separation in Figure 11. Sample names are given for some Type Ⅱ springs.

In addition, all deep groundwater samples were in equilibrium or supersaturated with calcite and dolomite, while some of the shallow groundwater was undersaturated with both (Figure 6A). Most of the Ca²⁺ seemed to originate from calcite dissolution in the shallow groundwater, while from dolomite dissolution in the deep groundwater (Figure 6B), according to Lacelle et al. (2008) who assumed that all Mg²⁺ in water is derived from dolomite dissolution and the relative contribution of Ca²⁺ from calcite dissolution can be estimated by subtracting the Mg²⁺ from Ca²⁺ concentrations in water.

FIGURE 6

FIGURE 6. Degrees of water-rock interactions: (A) Saturation index (SI) of calcite and dolomite; (B) Ca²⁺ derived from dolomite (estimated using mMg) vs. Ca²⁺ derived from calcite (estimated using mCa—mMg) following Lacelle et al. (2008). In (A) and (B), the slope of the solid line is one. (C) Ca²⁺/Mg²⁺ vs. Ca²⁺ + Mg²⁺; (D) Ca²⁺/(Na⁺ + K⁺) vs. Ca²⁺ + Mg²⁺ depending on the water group; (E) Ca²⁺/Mg²⁺ vs. Ca²⁺ + Mg²⁺; (F) Ca²⁺/(Na⁺ + K⁺) vs. Ca²⁺ + Mg²⁺ depending on the lithology around the sampling location. Rainwater was excluded in (E) and (F). In (D), the letters indicate the names of distinct samples.

Furthermore, the shallow groundwater had higher Ca²⁺/Mg²⁺ but lower Ca²⁺/(Na⁺ + K⁺) (Figures 6C,D) than the deep groundwater, probably due to the geology, as seen in Figures 6E,F, since most of the deep groundwater samples were obtained from limestone and dolomitic limestone and thus had high Ca²⁺ + Mg²⁺ and Ca²⁺/(Na⁺ + K⁺), whereas the shallow groundwater was obtained from various geological sources including the Bansong Formation (Figure 1; Supplementary Table S1). As the water-silicate interaction proceeded in the Bansong Formation, Ca²⁺/(Na⁺ + K⁺) decreased despite high Ca²⁺ + Mg²⁺ as Y37, MW3, and MW7 in Figures 6D,F.

Meanwhile, the surface water occurred geologically in the Yeongheung Formation (Oy in Figure 1) and had relatively high Mg²⁺ concentrations (Figure 4A), implying the hydraulic connectivity between surface water and groundwater as shown in the tracer test with uranine at the west side of the study area (out of the map in Figure 1A) in January 2021, and was saturated with calcite and dolomite (Figure 6A) except for a sample (ST4) from the Bansong Formation (Jbc in Figure 1), probably due to CO₂ degassing given low P_CO2 (average = 10^−2.93 atm in Table 1). A surface water sample (ST2) had a slightly low δ¹³C_DIC (−14.2‰) compared to the groundwater in Figure 5B, but still high compared to the biogenic CO₂ (−26 to −21‰ in Jiang et al., 2013; Kang et al., 2020), indicating the effect of other carbon sources (e.g., calcite).

Spring Water

Geographically, Type I springs mainly occurred in the Jurassic Bansong Formation occupying topographically high areas (220–610 m asl) and had an average water level of 333.5 m asl (n = 11), while Type II mainly discharged in the Yeongheung Formation and occasionally occurred in the Bansong Formation (Figures 1, 7) and had an average water level of 301.2 m asl (n = 30). Type Ⅱ was also hydrochemically distinct from Type Ⅰ by high EC, Mg²⁺, HCO₃⁻, and to a lesser degree, Ca²⁺ (Figure 4; Table 1). In addition, most of the Type Ⅱ springs were close to saturation with both calcite and dolomite (Figure 6A) and had high Ca²⁺ + Mg²⁺ and Ca²⁺/(Na⁺ + K⁺) (Figure 6D), while Type Ⅰ was undersaturated with both carbonate minerals (Figure 6A) and had most of the Ca²⁺ from calcite dissolution (Figure 6B) and low Ca²⁺ + Mg²⁺ and Ca²⁺/(Na⁺ + K⁺) (Figure 6D).

FIGURE 7

FIGURE 7. Major elements by lithology: Oy. = Yeongheung Formation (limestone and dolomitic limestone), Omg = Maggol limestone, Jbc = Bansong Formation (mainly, metamorphosed shale and sandstone) including Jbs in Figure 1A, Ogl = Kosong limestone (the white area at the east in Figure 1A). The letters on the diagram indicate the names of distinct samples.

The hydrochemical differences between Types Ⅰ and Ⅱ were conflicting with the expectation that conduit flow (Type Ⅱ in this study) is characterized by low contents of dissolved elements and is undersaturated with respect to calcite (e.g., De Rooij and Graham, 2017), which was probably due to the different geology (Figures 6F, 7) and water sources (Figure 8) as in Lorette et al. (2018) and supported by the PCA result (Figure 9A). Four principal components (PCs) were extracted, explaining 71% of the total variance (Table 2). PC1, which explained the highest fraction of explained variance, was highly correlated with EC, Ca²⁺, Mg²⁺, HCO₃⁻, SI_calcite, and SI_dolomite (Table 2), reflecting the dissolution of carbonate rocks. Most of the Type Ⅱ springs were plotted on the positive side of PC1 with deep groundwater (Figure 9A), implying the discharge of deep groundwater at Type Ⅱ springs in carbonate rocks (Figure 8), while Type Ⅰ was on the left side (Figure 9A), implying relatively little water-carbonate interaction (Figure 8). In addition, Type Ⅱ had relatively higher δ¹⁸O and δD values than Type Ⅰ in the dry season (Figure 5C) except for one sample (Y25-1 from the Kosong limestone in Figure 1A), and relatively higher temperature and lower DO (Table 1), which also supported the discharge of deep groundwater at Type Ⅱ (Figure 8).

FIGURE 8

FIGURE 8. Summary of subsurface flow in the studied karst area with Cambro-Ordovician dolomitic limestone overlying Jurassic metamorphosed shale and sandstone along thrusts. The conduit-flow-dominated Type Ⅱ springs discharge deep groundwater and contaminants infiltrating the subsurface, whereas the diffuse-flow-dominated Type Ⅰ springs discharge rainwater carrying contaminants and dissolving calcite via shallow circulating pathways and little water-rock interactions in mountain areas.

FIGURE 9

FIGURE 9. Principal component analysis: (A) scores of PC1 vs. PC2; (B) scores of PC1 vs. PC3. The letters on the diagram indicate the names of distinct samples.

TABLE 2

TABLE 2. Loadings of principal components (PCs). Absolute values of loadings ≥0.6 are in bold.

Similarly, Lorette et al. (2018) found high Mg²⁺ up to 12.8 mg/L and equilibrium with calcite during low-water periods at a Vauclusian-type spring located on a major faulted anticlinal structure, and explained the unusual hydrochemistry using the mixture with a Jurassic confined aquifer since high Mg²⁺ concentrations were only observed in Jurassic confined aquifers (> 20 mg/L as a result of the presence of magnesium-calcite) as in the study area (Table 1; Figure 3). Moore et al. (2009) suggested the deep-water source (Ca–Mg–SO₄-type water from Well 2 with Mg²⁺ up to 49.6 mg/L) and local diffuse recharge (Ca–HCO₃-type water from Well 4 with low Mg²⁺ between 1.2 and 2.2 mg/L resulting from rain water equilibrating with the Ocala limestone) in a system dominated by conduit flow to explain the changes in chemical compositions (e.g., Na⁺, Mg²⁺, K⁺, Cl⁻, and SO₄²⁻) of spring discharge. Low concentrations of dissolved ions (e.g., Mg²⁺ and SO₄²⁻) were explained with little contribution of deep water (Moore et al., 2009). Bicalho et al. (2012) concluded high and low mineralized groundwaters with Cl⁻, Na⁺, Mg²⁺, and SO₄² respectively, related to deep groundwater rising and to superficial infiltration water at the Lez spring in the Mediterranean basin, and showed a decrease of the deep compartment participation to the Lez spring outflow due to intense exploitation that decreases the total hydraulic head in the aquifer. Fiorillo et al. (2018) deducted upwelling groundwater (with high ²²²Rn activity) flux feeding the springs by the increasing of the hydraulic head in depth. The high flow condition is sustained by high ascendant hydraulic gradient, under which higher water velocity and higher radon activity occur.

Meanwhile, spring water contained high NO₃⁻ compared to groundwater (Figures 3, 4C). The contamination of spring water with NO₃⁻ was supported by PC3 which was positively related with Cl⁻ and NO₃⁻ (Table 2), indicating anthropogenic contamination given that Cl⁻ and NO₃⁻ are major contaminants in agricultural areas (e.g., Schilling and Helmers, 2008; Kim et al., 2021). Most of the groundwater was plotted on the negative side of PC3, and a shallow groundwater well (Y37-1) occurring on the Bansong Formation with 0.2–0.8 mg/L for NO₃⁻ and 1.6–2.1 mg/L for Cl⁻ had the lowest PC3 scores (Figure 9B). In contrast, about half of the spring water samples were on the positive side of PC3. In particular, Type Ⅱ had high NO₃⁻ and Cl⁻ concentrations (Figure 4C), implying the fast infiltration of contaminants through conduits (Figure 8). The fast infiltration from surface was supported by the tritium concentrations (Figure 5A) similar to the ³H levels around 10 TU in precipitation during the study period in South Korea (Chae and Kim, 2019). Type Ⅰ also showed relatively high tritium concentrations but less than Type Ⅱ, which suggests the recent or short-circulating discharge of meteoric water occurring at both springs (Palcsu et al., 2021).

In summary, the Type Ⅱ springs seemed to be fed by both deep groundwater and recent conduit flow (Figure 8) as in Lorette et al. (2018) who suggested that Toulon Springs are fed by two aquifers: (i) the Jurassic confined karst aquifer responsible for low-water period and related to dolomitic hydrochemical characteristics, and (ii) the upper Cretaceous unconfined aquifer responsible for quick hydrodynamic and hydrochemical variations. In contrast, Type Ⅰ seemed to discharge diffuse sources with low residence time (Figure 8) as in the Tufa Springs system provided by Tobin and Schwartz (2012) who calculated specific conductance in unsampled diffuse sources in mountain hydrologic systems that ranged from 34 to 257 mS/cm.

δ¹³C_DIC was in a wide range (−15.3–−11.2‰ in Figure 5B), implying multiple carbon sources from calcite, atmospheric CO₂, and soil CO₂ given the δ¹³C of carbonate rocks (∼0.4‰), atmospheric CO₂ (−9–−7‰), and soil CO₂ (−26–−21‰) (Jiang et al., 2013; Kang et al., 2020), and land use (Supplementary Figures S1, S2).

Major Hydrochemical Processes

Based on the PCA results (Table 2 and Figure 9) and hydrochemical and isotopic characteristics, water-rock interactions with carbonate (PC1) and silicate rocks (PC2) and anthropogenic contamination (PC3) can be suggested as major hydrochemical processes in the study area. PC4 was positive with pH but negative with SiO₂, and surface water and some deep groundwater samples had high PC4 scores, whereas PC4 could not be clearly explained given the information and thus was excluded from further discussion.

Effects of Geology

Type Ⅱ springs and deep groundwater had high EC, Mg²⁺, HCO₃⁻, SI_calcite, and SI_dolomite, which were positively correlated with PC1 (Table 2), and also had high Ca²⁺ + Mg²⁺ and Ca²⁺/(Na⁺ + K⁺) (Figure 6D). These properties can be explained by the carbonate dissolution given that most of Type Ⅱ springs and deep groundwater occurred in the limestone and dolomitic limestone (Figures 1, 6F, 7), whereas a surface water sample (ST4) and Type Ⅰ spring sample (S68-3) with negative PC1 scores occurred in the Bansong Formation (Figure 1) and were undersaturated with calcite and dolomite (Figure 6A).

However, the deep groundwater had Ca²⁺ from dolomite dissolution (Figure 6B) and low Ca²⁺/Mg²⁺ close to one (Figure 6C), whereas the Type Ⅱ springs contained higher Ca²⁺ than the deep groundwater (Table 1; Figure 4B), probably because of the extensive water-rock interaction in deep groundwater, causing calcite precipitation (Figure 6A) and thus decreasing Ca²⁺/Mg²⁺. According to Langmuir (1971), the molar Ca²⁺/Mg²⁺ ratio as low as 0.6 results from incongruent dissolution of dolomite. Besides, calcite is precipitated at low temperature when the molar Ca²⁺/Mg²⁺ ratio is greater than one given a large uncertainty (Drever, 1997). The deep groundwater showed the molar Ca²⁺/Mg²⁺ ratios between 0.8 and 1.4 (average of 1.1; Figure 6C), which indicates the possible condition for the incongruent dissolution of dolomite and calcite precipitation. The lower average P_CO2 of deep groundwater (10^−2.46 atm) than that of Type Ⅱ springs (10^−2.28 atm) also supported the prolonged dissolution of calcite in deep groundwater in a closed system (Figure 10) since CO₂ is consumed for calcite dissolution (Langmuir, 1971; Lacelle et al., 2008; Frondini et al., 2019). Meanwhile, the inconsistent Mg²⁺ contents in carbonate rocks might affect the Ca²⁺/Mg²⁺ ratios. For instance, the water sample from the Kosong limestone had high Ca²⁺/Mg²⁺ regardless of water groups (Figure 6E), probably because the Kosong limestone consists mainly of calcite.

FIGURE 10

FIGURE 10. (A) Saturation index (SI) of calcite vs. P_CO2. (B) A diagram of pH versus HCO₃⁻ showing evolution as water dissolves calcite under open (solid) and closed (dotted) system conditions with varying initial P_CO2.

PC2 was positively correlated with Na⁺ and K⁺ but negatively with Eh (Table 2). Some groundwater samples had high positive PC2 scores, while most springs were plotted on the negative side of PC2 (Figure 9A), indicating that PC2 addressed silicate dissolution. In particular, a deep groundwater sample (MW7) from the Bansong Formation had the highest PC2 score. This deep groundwater sample had the highest Na⁺ and F⁻ (Figure 7). In addition, some groundwater had relatively low Ca²⁺/(Na⁺ + K⁺) values due to high Na⁺ and K⁺ in the Bansong Formation (e.g., Y37, MW3, and MW7; Figures 6F, 7), indicating silicate dissolution and thus geological heterogeneity in this carbonate aquifer. Especially, two deep groundwater samples (MW7, MW3) imply extensive water-silicate interaction but little anthropogenic contamination in the deep aquifer based on low Ca²⁺/(Na⁺ + K⁺) (Figure 6D), high Na⁺ and F⁻ concentrations (Figure 7), and low NO₃⁻ and Cl⁻ concentrations (Figure 4C).

Note that Type Ⅰ springs also had low Ca²⁺/(Na⁺ + K⁺) values and occurred in the Bansong Formation (Figures 6D,F). However, the effect of silicate dissolution seemed slight in Type Ⅰ given the similar concentrations of Na⁺, K⁺, and F⁻ in both types of springs (Table 1) despite the different geology (Figures 1, 7). The low contents of Na⁺ and K⁺ as well as EC and Ca²⁺ would be probably because of the discharge of diffuse sources with low residence time at Type Ⅰ springs (Figure 8) as in the mountain systems provided by Tobin and Schwartz (2012), which also can explain Type Ⅰ containing lower concentrations than the shallow groundwater for all major ions except NO₃ (Table 1). Similarly, Moore et al. (2009) showed low and relatively constant concentrations (e.g., Mg²⁺ and SO₄²⁻) of diffuse recharge at a well in a karstic aquifer due to low inputs of deep water.

Anthropogenic Contamination

Spring water contained higher NO₃ concentrations (Figure 4C) and PC3 scores (Figure 9B) than groundwater, indicating that the contaminants stored in the vadose zone discharged through springs. In particular, the point infiltration through conduits seemed to deliver more contaminants, causing higher Cl⁻ and NO₃⁻ concentrations at Type Ⅱ (Figure 4C) given that contaminant storage occurs in the rock matrix and epikarst, while contaminants are mostly transported along preferential pathways in karst systems (Ghasemizadeh et al., 2012). It should be noted that the Type Ⅰ springs occasionally (e.g., S37) had higher ratios of Cl⁻ + NO₃⁻ than Type Ⅱ on the Piper’s diagram in Figure 3, which was because the Type Ⅰ springs contained relatively low HCO₃⁻ concentrations similar to rainwater (Figure 4A) implying the discharge of diffuse sources with little water-rock interaction (Figure 8).

Shallow groundwater was also affected by the influx of anthropogenic contaminants from the surface, given the high NO₃⁻ or Cl⁻ concentrations (Figure 4C). In particular, the Y28 had the highest PC3 score (Figure 9B) due to a high level of both NO₃⁻ and Cl⁻ (Figure 4C) and showed a high level of Mg²⁺, Na⁺, and SO₄²⁻ (Figure 7), which indicates the combined effect of carbonate and silicate dissolution and anthropogenic contamination since the subsurface can be the Bansong Formation at the geologic boundary, although the surface geology of Y28 was the Yeongheung Formation (Figure 1). Also, the sampling site (Y28) was close to the Dong River and agricultural land (Supplementary Figure S2). In addition, Y37 had high K⁺ and Cl⁻ (Figure 7) and thus low Ca²⁺/(Na⁺ + K⁺) (Figure 6D). MW2 also had high K⁺ (Figure 7). High K⁺, Na⁺, and F⁻ were occasionally found in groundwater (Figure 7) and the average Na⁺ and K⁺ concentrations of groundwater were larger than the corresponding concentration of spring water (Table 1), implying the effect of silicate dissolution based on their (e.g., MW7, Y37) occurrence in the Bansong Formation. However, the effect of agrochemicals cannot be ignored for high K⁺ as in Stueber and Criss (2005) and Moore et al. (2009) based on land use (Supplementary Figure S2) and the fact that K⁺ increased with Cl⁻ during the wet season (July) in shallow groundwater (e.g., Y37 and MW2 in Figure 7 and Supplementary Figure S3).

Furthermore, HCO₃⁻ concentrations were higher than those predicted from calcite dissolution at a given P_CO2 (Supplementary Figure S4a), which indicates additional proton (H⁺) sources. There was no evidence of pyrite to produce H⁺ in the study area, while nitrate contamination might enhance carbonate dissolution (Kim et al., 2019). In fact, the water samples were plotted between the enhanced dissolution of carbonate rocks by nitrification (i.e., the 1:1 of Ca²⁺ + Mg²⁺ (mmol) and HCO₃⁻ (mmol)) and natural dissolution (i.e., the 1:2 of Ca²⁺ + Mg²⁺ (mmol) and HCO₃⁻ (mmol)) (Supplementary Figure S4b). PC3 was negatively related with pH but positively with NO₃⁻ (Table 2), implying the acidification and nitrate contamination of water.

Subsurface Flow in the Study Area

Based on the hydrochemical difference of spring water depending on the flow type and PC1 (carbonate dissolution) in Figure 9A and their spatial distribution on the geologic map in Figure 1A, it can be deduced that the heterogeneous geology causes different flow patterns in the study area. The carbonate rocks seemed to provide the conduits for Type Ⅱ due to their relatively high solubility and fractures (Figure 1B), whereas the sandstone and shale provided a matrix for the diffusive flow of Type Ⅰ (Figure 8). The different hydrological conditions may further differentiate the hydrochemistry, by causing different contributions of end-members (e.g., groundwater, rainfall) in spring discharge, as in Bicalho et al. (2012) who found that hydrological conditions induce the different proportion of end-members (deep groundwater rising and superficial infiltration water) participating in groundwater mixing of a spring.

Type Ⅰ Spring

Type Ⅰ springs seemed to slowly discharge the rainwater infiltrated from soil cover in the Bansong Formation, carrying contaminants (e.g., NO₃ and Cl) and dissolving calcite via shallow (near-surface) circulating pathways (Figure 8) based on the following results: low concentrations of major compositions including HCO₃⁻ and EC (Figure 4; Table 1), undersaturation with respect to calcite (Figure 6A) and Ca²⁺ from calcite dissolution (Figure 6B), and little water-rock interaction (i.e., negative PC1 and PC2 scores in Figure 9A). Low P_CO2 similar to that in rainwater despite slow flow also suggests a short residence time or bare soil given soil P_CO2 around 10^−1.5–10^−2.5 atm (Appelo and Postma, 2005) or higher due to an accumulation of CO₂ in the epikarst (Peyraube et al., 2012). Type Ⅰ seemed to represent springs in mountains following a relatively shallow flow path with little contribution of groundwater because of steep slopes and shallow soil development (Somers and McKenzie, 2020).

Most of Type Ⅰ springs were undersaturated with respect to calcite and showed P_CO2 near typical biologically respired values (10^−2.5 atm in Lacelle et al., 2008; Figure 10A). In addition, the plotting of the Type Ⅰ springs on pH against HCO₃⁻ suggests an evolutionary path by the dissolution of carbonates under open system conditions with P_CO2 near 10^−2.5 atm (Figure 10B), which increases both the pH and HCO₃⁻ of spring water with subsurface residence times. The significant amounts of NO₃⁻ and high DO also suggested an open system (Table 1).

Type Ⅱ Spring

The relatively high levels of NO₃⁻, Cl⁻ (Figure 4C), tritium, and low δ¹³C_DIC down to –15.3‰ in Type Ⅱ springs (Figure 5) supported the fast flow of recent or short-circulating recharge through conduits and the admixture of biologically sourced carbon (Figure 8). In addition, Type Ⅱ had high Mg²⁺, HCO₃⁻ (Figure 4A), δ¹⁸O, δD (Figure 5C), SI_calcite, and SI_dolomite (Figure 6A) similar to the deep groundwater, which can be explained by the discharge of deep groundwater at Type Ⅱ springs (Figure 8). Other possible causes such as CO₂ exsolution, increases in temperature, or dissolution of soluble minerals were excluded given the higher P_CO2 (average = 10^−2.28 atm) and lower temperature (14.7°C) in Type Ⅱ springs than those in groundwater (Table 1) and the fast flow rates (Figure 2). A negative relation between P_CO2 and SI_calcite (Figure 10A) and the relationship between pH and HCO₃⁻ along the saturation curve (Figure 10B) also support the hydrochemical changes of Type Ⅱ springs due to changes in water chemistry that affect P_CO2 (e.g., carbonate dissolution or precipitation), similar to the deep groundwater.

Type Ⅱ seemed to discharge both the flow component through the epikarst zone and the baseflow from the lower reservoir, as defined by Tritz et al. (2011) who used two reservoirs (i.e., the upper epikarst/soil zone and the lower vadose and saturated zone) to simulate the behavior of a karst system catchment and assumed a fast flow through the epikarst and slow baseflow from the lower reservoir. The deep groundwater rising at springs was also observed by other researchers including Moore et al. (2009), Bicalho et al. (2012), Demiroglu (2016), Gil-Márquez et al. (2019), Lorette et al. (2018), and Zhu et al. (2020). Some implied fast transfer through conduits unlike Tritz et al. (2011). For instance, Lorette et al. (2018) showed Toulon Springs with a mean annual daily discharge of 450 L/s fed by deep confined karst aquifers during low-water periods. Gil-Márquez et al. (2019) mentioned the piston-flow effect at the beginning of floods, which causes the drainage of deep ascending flows through the saturated zone and subsequently the rise in mineralization at a spring. Baudement et al. (2017) mentioned chimney-shafts, by which deep phreatic karst systems are connected to vauclusian springs. Based on these previous study results and the results of this study, it can be inferred that karst conduit networks including palaeo-karstic systems or abandoned cavities seemed to still be communicating with the main drainage axis in the study area as in Aquilina et al. (2005) and Pogačnik et al. (2017).

Hydrograph Separation of Type Ⅱ Springs

Proportion of Groundwater at Type Ⅱ Springs

The proportion of baseflow at Type Ⅱ springs was evaluated using the hydrograph separation. δ¹⁸O and δD decreased as it began to rain heavily at three Type Ⅱ springs and increased back (Figure 11A), reflecting the changing of water sources (e.g., rainwater (R) and pre-storm water (PS)) in the spring discharge during the typhoon event. The PS was represented using the average δD and δ¹⁸O of deep groundwater collected during the dry season (δ¹⁸O = −9.6 and δD = −70.2‰; n = 3; Figure 5C), given the hydrochemistry of Type Ⅱ similar to that of deep groundwater. It was assumed that the deep groundwater in the dry season records the isotopic compositions of baseflow at Type Ⅱ. In fact, the spring water collected during the same period (n = 11; Figure 5C) had a similar isotopic composition to the groundwater on average, indicating baseflow discharge in the dry season. All the δD and δ¹⁸O were plotted slightly below the meteoric water line (Figure 5C), indicating evaporation or water mass mixing (Serno et al., 2017). Also, all three springs had initial values similar to the average of deep groundwater in Figures 5D, 11A, reflecting baseflow discharge before the storm.

FIGURE 11

FIGURE 11. Isotropic variations (A) and hydrograph separation (B) at three Type Ⅱ springs during the typhoon in Figure 2B. The blue bar indicates precipitation.

In contrast, a rainwater sample (R) collected during the storm event had −12.0‰ for δ¹⁸O and −85.0‰ for δD (Figure 5D). These values were quite different from those of PS. This large isotopic difference in the two end-members (i.e., 2.4‰ for δ¹⁸O and 14.3‰ for δD) and the plotting of spring water between R and PS supported the applicability of the two-component hydrograph separation in the study area. Given the low δ¹⁸O and δD of R (Figure 5D), the lower δ¹⁸O and δD value of spring water indicates a higher proportion of R in discharge; thus, high proportions of R were expected in discharge at S46 given the great decreases in δ¹⁸O and δD (Figure 11A).

We chose the bigger proportion of R between those using δ¹⁸O and δD (Table 3; Figure 11B). As a result, the proportion of R was similarly 17, 13, and 13% of the total discharge at S41, S45, and S46, respectively, before the storm. At peak flow, however, R accounted for 53, 39, and 87% at S41, S45, and S46, respectively, while the maximum contribution of R was 60, 53, and 87%.

TABLE 3

TABLE 3. Contribution (%) of rainfall (R) to discharge at three Type Ⅱ springs during the storm event in Figure 2B. The remaining was assumed to come from pre-storm water (PS; deep groundwater). The bigger proportion of R between those using δ¹⁸O and δD was in bold. Time for peak flow was shaded in gray.

This result indicates that the proportion of groundwater was up to 83–87% in the three Type Ⅱ springs in the dry season, similar to Moore et al. (2009) who identified that upwelling from deep flow paths provides significant contributions of water to spring discharge and emphasized that restricted monitoring of springs limits the interpretation of karst systems by masking critical components of the aquifer. The proportion of baseflow decreased during the storm, while the contribution of rainwater increased at Type Ⅱ (Table 3; Figure 11) as in Moore et al. (2009) who showed that the contribution from the deep source depends inversely on flow conditions.

Note the high tritium concentrations at Type Ⅱ (Figure 5A) despite the high proportion of deep groundwater in the dry season was probably because the water samples for tritium were obtained in wet season including S46.

Storage Effect

S46 had the maximum contribution of rainfall (87%) at peak flow without a delay, while S41 and S45 had a delay between the peak flow and maximum contribution of R (Table 3). At the end of monitoring, the proportion of R decreased to 21 and 17% at S41 and S45, respectively, indicating the recover to the baseflow condition, whereas it was still 60% at S46 unlike before the storm (i.e., 13%).

In addition, S46 did not show an increase in NO₃⁻ during the typhoon in 2000. Rather, the NO₃⁻ was diluted as other major compositions during the wet season in 1999 (Figure 12; Supplementary Figure S5). In contrast, NO₃⁻ increased at S41 and S45 during the typhoon event in 2000 as well as wet season in 1999 compared to during the dry season in 1999, indicating the flushing out of NO₃ stored in the vadose zone, as observed by other researchers including Doctor et al. (2006) and Lorette et al. (2018). The influx of NO₃⁻ from the vadose zone seemed to be more influential than the dilution effect of rainwater at S41 and S45, as in Lorette et al. (2018) who observed that during recharge events the first flood event leads to an increase of NO₃⁻ concentrations, resulting in the mobilization of NO₃⁻ storage in soil and the unsaturated zone of the aquifer.

FIGURE 12

FIGURE 12. The temporal variation of major compositions at three Type Ⅱ springs in the precipitation events in Figure 2. The data in 2000 were the average of 10 samples collected for 5 days (Kim et al., 2001).

Given the similar NO₃ levels (Figure 12) as well as geographical locations (Figure 1), geology (Yeongheung Formation), and flow type (Type Ⅱ), the distinguished behavior of NO₃⁻ (Figure 12) and water isotopes (Figure 11) at S46 was probably due to the low storage capacity of the vadose zone for water and contaminants. Water and solutes stored in the vadose zones can be discharged through springs following rainfall events in karst areas (Perrin et al., 2003; Aquilina et al., 2005), which may cause buffered rain isotope responses despite highly peaking hydraulic response (Perrin et al., 2003), increase P_CO2 and anthropogenic components (e.g., Cl⁻) and decrease δ¹³C_DIC by the contribution of water held in storage in the vadose zone to the conduit system in storm events (Doctor et al., 2006), or increase NO₃⁻ at the first flood event during recharge events, with the decrease in NO₃⁻ in next flood events due to the stock of NO₃⁻ storage in the soil and the unsaturated zone moved out of the karst system during the previous flood event (Lorette et al., 2018). Jacob et al. (2008) showed spatial and temporal variations of vadose zone water storage in a karst system using ground-based absolute gravimetry measurements.

The low storage seemed to cause the discharge quickly responding to rainfall not only during the typhoon event in 2000, but also in the small precipitation events in 1999 (Figure 2), and the rainfall accounting for 87% of the discharge at peak flow and continuing to flow out and contributing more than half of the discharge after the heavy rainfall stopped at S46 (Figure 11; Table 3). Similarly, Demiroglu (2016) reported that developed karst sinkholes allow fast percolation of heavy rainfall into an aquifer up to 80%, and explained that the very low storage combined with the high transmissivity means that most of the recharge will not be retained by the karst system, but will rapidly flow out to springs, rivers, lakes, or seas. In contrast, the delayed response to the storm and quick recovery to the baseflow condition after it stopped raining at S41 and S45 (Figure 11; Table 3) were probably due to the storage effect of the vadose zone, given that all three springs had flow characteristics through permeable conduits. Doctor et al. (2006) considered the anthropogenic component stored within the epikarst as one of the recharge sources, and found that an anthropogenic component derived from epikarstic storage affects the well under conditions of elevated hydraulic head, accounting for the chemical response in the well during wet conditions. A similar increase in NO₃⁻ and Cl⁻ during the wet season was observed at the wells L6-1 (deep) and Y37-1 (shallow) as well as the increases in K⁺ and Cl⁻ at MW2 and Y37 (shallow) in the study area (Supplementary Figure S3).

We acknowledged that the storage effect of the subsurface in the study area needs to be further investigated since the storage degree was evaluated only with the response process at this current study, although the subsurface storage is involved with many processes, including response, recession, and discharge, which are also affected by aquifer structures such as the development of conduits, degree of weathering, and recharge areas. In addition, there is an opposite interpretation from Demiroglu (2016) who mentioned that the variability of physicochemical characteristics is higher in conduit permeability with low storage. Stueber and Criss (2005) suggested that less soluble ions increase with flow as they are mobilized from fields to karst conduits under storm conditions, while highly soluble ions supplied by diffuse groundwater are diluted by high flows.

Conclusion

Spring water quality was mainly controlled by three hydrogeological processes in the studied karst area: 1) carbonate dissolution (PC1 in Figure 9) anthropogenic contamination (PC3 in Figure 9) rainfall infiltration (Figure 11). Precisely, the Type Ⅱ springs, which was assumed to be conduit flow-dominated, were greatly affected by carbonate dissolution and had high EC, Ca²⁺, Mg²⁺, HCO₃⁻, Ca²⁺ + Mg²⁺, Ca²⁺/(Na⁺ + K⁺), and SI_calcite, similar to the deep groundwater, probably because of their spatial occurrence in the carbonate rocks and Type Ⅱ springs supplied by the deep groundwater (Figure 8). Hydrograph separation using δ¹⁸O and δD data showed that the deep groundwater contributed to 83–87% of the discharge in three Type Ⅱ springs in the dry season. The contribution of baseflow decreased as the contribution of rainfall increased during storm events at Type Ⅱ springs. A Type II spring (S46) rapidly responded to rainfall and showed the contribution of rainwater exceeding 50% 2.5 days after the rain stopped, with little NO₃ released from the vadose zone during the storm, implying low storage. In contrast, the other two Type II springs (S41 and S45) had delayed responses to precipitation and had NO₃ flushed out from the vadose zone during the wet season, indicating that the vadose zone may store rainwater and release contaminants stored depending on the type of precipitation.

In contrast, the effect of silicate dissolution (PC2 in Figure 9) appeared in groundwater but not in spring water, probably due to the short circulation of Type I springs despite their occurrence in the metamorphosed sandstone and shale. Type Ⅰ springs seemed to slowly discharge infiltrated rainwater, dissolving carbonate minerals in the relatively high-altitude area. Both types of springs were contaminated with NO₃⁻ and Cl⁻ compared to the groundwater, and in particular, Type Ⅱ had high average concentrations of Cl⁻ and NO₃⁻, indicating severe contamination through conduits.

Based on the study results, the heterogeneous geology affected both the hydrochemistry and flow types of spring water in this stratigraphically complex karst area with thrusts, causing unusual hydrochemical characteristics: elevated concentrations of Mg²⁺ at Type Ⅱ and low concentrations at Type Ⅰ. The carbonate rocks provided the conduits for Type Ⅱ, whereas the silicate rocks provided a matrix for Type Ⅰ, which further caused hydrochemical differences depending on the flow type.

The study results suggest that hydrochemical and isotopic data are applicable to understand the flow paths of springs (e.g., conduit flow through carbonate rocks) and the storage characteristics in the vadose zone in a karst area since the hydrological condition (e.g., upwelling of deep groundwater, rainfall recharge through a matrix) affects the participation of water sources and consequently the hydrochemical and isotopic composition of springs. Mg²⁺, Ca²⁺/Mg²⁺, Ca²⁺/(Na⁺ + K⁺), and SI_calcite were good indicators to distinguish the geology that the water passed through (e.g., limestone, dolomitic limestone, sandstone, and shale) and the hydrochemical evolution of water samples (e.g., incongruent dissolution of dolomite) in the tectonically complex karst area. PCA was useful to characterize the hydrochemical processes in each type of spring. Also, the present study showed the successful applicability of water isotopic tracers for evaluating the hydrodynamic properties of subsurface groundwater flow, although this study failed to quantify the effect of epikarst or soil water compared to rainfall and groundwater during storm events due to little information, which remains future work. Besides, upwelling of deep groundwater at Type Ⅱ springs estimated by the hydrochemical and isotopic analysis should be confirmed by a hydrogeological investigation. In particular, the main drainage network for deep groundwater upwelling needs to be identified in the near future given a new plan to build a waste landfill in a historical limestone mine site around the study area. Thrust faults should be first assessed as a potential pathway. In addition, the classification of springs is needed to be confirmed based on Reynolds numbers.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

Author Contributions

SY interpreted the data and wrote the manuscript, GC collected the data and drew figures, JO conducted the principal component analysis and drew figures, S-HK and D-IK collected and interpreted the data, and S-TY initiated the project and wrote the manuscript.

Funding

The initiation of this work was funded by the Hyundai Engineering Co. as a part of a site investigation for a dam in 1999–2000. The completion of this work was supported by the Korea Ministry of Environment (MOE) as the “Korea-CO₂ Storage Environmental Management (K-COSEM) Research Program.” The first author was also supported by a National Research Foundation of Korea (NRF) grant funded by the Korea government (MEST) (No. 2019R1A2C1084297).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We would like to thank Jieun Seo for advice on the geology of the study area.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/feart.2021.712865/full#supplementary-material

References

Adji, T. N., Haryono, E., Fatchurohman, H., and Oktama, R. (2016). Diffuse Flow Characteristics and Their Relation to Hydrochemistry Conditions in the Petoyan Spring, Gunungsewu Karst, Java, Indonesia. Geosci. J. 20 (3), 381–390. doi:10.1007/s12303-015-0048-8

CrossRef Full Text | Google Scholar

Appelo, C. A. J., and Postma, D. (2005). Geochemistry, Groundwater and Pollution. 2nd Edition. CRC Press.

Aquilina, L., Ladouche, B., and Dörfliger, N. (2005). Recharge Processes in Karstic Systems Investigated Through the Correlation of Chemical and Isotopic Composition of Rain and Spring-Waters. Appl. Geochem. 20, 2189–2206. doi:10.1016/j.apgeochem.2005.07.011

CrossRef Full Text | Google Scholar

Aquilina, L., Ladouche, B., and Dörfliger, N. (2006). Water Storage and Transfer in the Epikarst of Karstic Systems During High Flow Periods. J. Hydrol. 327, 472–485. doi:10.1016/j.jhydrol.2005.11.054

CrossRef Full Text | Google Scholar

Atkinson, T. C. (1977). Diffuse Flow and Conduit Flow in Limestone Terrain in the Mendip Hills, Somerset (Great Britain). J. Hydrol. 35, 93–110. doi:10.1016/0022-1694(77)90079-8

CrossRef Full Text | Google Scholar

Barbieri, M., Boschetti, T., Petitta, M., and Tallini, M. (2005). Stable Isotope (2H, 18O and 87Sr/86Sr) and Hydrochemistry Monitoring for Groundwater Hydrodynamics Analysis in a Karst Aquifer (Gran Sasso, Central Italy). Appl. Geochem. 20, 2063–2081. doi:10.1016/j.apgeochem.2005.07.008

CrossRef Full Text | Google Scholar

Baudement, C., Arfib, B., Mazzilli, N., Jouves, J., Lamarque, T., and Guglielmi, Y. (2017). Groundwater Management of a Highly Dynamic Karst by Assessing Baseflow and Quickflow With a Rainfall-Discharge Model (Dardennes Springs, SE France). Bull. Soc. Géol. Fr. 188, 40. doi:10.1051/bsgf/2017203

CrossRef Full Text | Google Scholar

Bicalho, C. C., Batiot-Guilhe, C., Seidel, J. L., Van Exter, S., and Jourde, H. (2012). Hydrodynamical Changes and Their Consequences on Groundwater Hydrochemistry Induced by Three Decades of Intense Exploitation in a Mediterranean Karst System. Environ. Earth Sci. 65, 2311–2319. doi:10.1007/s12665-011-1384-2

CrossRef Full Text | Google Scholar

Chae, J.-S., and Kim, G. (2019). Seasonal and Spatial Variations of Tritium in Precipitation in Northeast Asia (Korea) Over the Last 20 Years. J. Hydrol. 574, 794–800. doi:10.1016/j.jhydrol.2019.04.058

CrossRef Full Text | Google Scholar

Chang, Y., Wu, J., Jiang, G., Liu, L., Reimann, T., and Sauter, M. (2019). Modelling Spring Discharge and Solute Transport in Conduits by Coupling CFPv2 to an Epikarst Reservoir for a Karst Aquifer. J. Hydrol. 569, 587–599. doi:10.1016/j.jhydrol.2018.11.075

CrossRef Full Text | Google Scholar

Chough, S., Kwon, S.-T., Ree, J.-H., and Choi, D. K. (2000). Tectonic and Sedimentary Evolution of the Korean Peninsula: a Review and New View. Earth-sci. Rev. 52, 175–235. doi:10.1016/s0012-8252(00)00029-5

CrossRef Full Text | Google Scholar

De Filippi, F., Iacurto, S., Grelle, G., and Sappa, G. (2021). Magnesium as Environmental Tracer for Karst Spring Baseflow/Overflow Assessment-A Case Study of the Pertuso Karst Spring (Latium Region, Italy). Water. 13, 93. doi:10.3390/w13010093

CrossRef Full Text | Google Scholar

De Rooij, R., and Graham, W. (2017). Generation of Complex Karstic Conduit Networks With a Hydrochemical Model. Water Resour. Res. 53, 6993–7011. doi:10.1002/2017WR020768

CrossRef Full Text | Google Scholar

Demiroglu, M. (2016). Classification of Karst Springs for Flash-Flood-Prone Areas in Western Turkey. Nat. Hazards Earth Syst. Sci. 16, 1473–1486. doi:10.5194/nhess-16-1473-2016

CrossRef Full Text | Google Scholar

Doctor, D. H., Alexander, E. C., Petrič, M., Kogovšek, J., Urbanc, J., Lojen, S., et al. (2006). Quantification of Karst Aquifer Discharge Components During Storm Events through End-Member Mixing Analysis Using Natural Chemistry and Stable Isotopes as Tracers. Hydrogeol. J. 14, 1171–1191. doi:10.1007/s10040-006-0031-6

CrossRef Full Text | Google Scholar

Drever, J. I. (1997). The Geochemistry of Natural Waters. 3rd edition. Englewood Cliffs, N.J: Prentice-Hall, 436.

Fiorillo, F., Esposito, L., Testa, G., Ciarcia, S., and Pagnozzi, M. (2018). The Upwelling Water Flux Feeding Springs: Hydrogeological and Hydraulic Features. Water. 10, 501. doi:10.3390/w10040501

CrossRef Full Text | Google Scholar

Frank, S., Goeppert, N., Ohmer, M., and Goldscheider, N. (2019). Sulfate Variations as a Natural Tracer for Conduit‐Matrix Interaction in a Complex Karst Aquifer. Hydrological Process. 33, 1292–1303. doi:10.1002/hyp.13400

CrossRef Full Text | Google Scholar

Frondini, F., Vaselli, O., and Vetuschi Zuccolini, M. (2019). Consumption of Atmospheric Carbon Dioxide through Weathering of Ultramafic Rocks in the Voltri Massif (Italy): Quantification of the Process and Global Implications. Geosciences. 9, 258. doi:10.3390/geosciences9060258

CrossRef Full Text | Google Scholar

Ghasemizadeh, R., Hellweger, F., Butscher, C., Padilla, I., Vesper, D., Field, M., et al. (2012). Review: Groundwater Flow and Transport Modeling of Karst Aquifers, With Particular Reference to the North Coast Limestone Aquifer System of Puerto Rico. Hydrogeol. J. 20, 1441–1461. doi:10.1007/s10040-012-0897-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Gil-Márquez, J. M., Andreo, B., and Mudarra, M. (2019). Combining Hydrodynamics, Hydrochemistry, and Environmental Isotopes to Understand the Hydrogeological Functioning of Evaporite-Karst Springs. An Example From Southern Spain. J. Hydrol. 576, 299–314. doi:10.1016/j.jhydrol.2019.06.055

CrossRef Full Text | Google Scholar

Grasso, D. A., Jeannin, P.-Y., and Zwahlen, F. (2003). A Deterministic Approach to the Coupled Analysis of Karst Springs' Hydrographs and Chemographs. J. Hydrol. 271, 65–76. doi:10.1016/s0022-1694(02)00321-9

CrossRef Full Text | Google Scholar

Han, R., Ree, J.-H., Cho, D.-L., Kwon, S.-T., and Armstrong, R. (2006). SHRIMP U-Pb Zircon Ages of Pyroclastic Rocks in the Bansong Group, Taebaeksan Basin, South Korea and Their Implication for the Mesozoic Tectonics. Gondwana Res. 9, 106–117. doi:10.1016/j.gr.2005.06.006

CrossRef Full Text | Google Scholar

Hartmann, A., Goldscheider, N., Wagener, T., Lange, J., and Weiler, M. (2014). Karst Water Resources in a Changing World: Review of Hydrological Modeling Approaches. Rev. Geophys. 52, 218–242. doi:10.1002/2013rg000443

CrossRef Full Text | Google Scholar

Jacob, T., Bayer, R., Chery, J., Jourde, H., Moigne, N. L., Boy, J.-P., et al. (2008). Absolute Gravity Monitoring of Water Storage Variation in a Karst Aquifer on the Larzac Plateau (Southern France). J. Hydrol. 359, 105–117. doi:10.1016/j.jhydrol.2008.06.020

CrossRef Full Text | Google Scholar

Jiang, Z., Lian, Y., and Qin, X. (2013). Carbon Cycle in the Epikarst Systems and its Ecological Effects in South China. Environ. Earth Sci. 68, 151–158. doi:10.1007/s12665-012-1724-x

CrossRef Full Text | Google Scholar

Kang, Y.-J., Yun, S.-T., Yu, S., Do, H.-K., and Chae, G. (2020). Quantitative Assessment of Deep-Seated CO2 Leakage Around CO2-Rich Springs with Low Soil CO2 Efflux Using End-Member Mixing Analysis and Carbon Isotopes. J. Environ. Manage. 276, 111333. doi:10.1016/j.jenvman.2020.111333

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, D. I., Yun, S. T., and So, C. S. (2001). “Impact of Rainfall on the Hydrochemical and Isotopic Characteristics in Spring Waters in a Karst Area,” in Proceedings of the Korean Society of Economic and Environmental Geology Conference, 23–26. (In Korean).

Google Scholar

Kim, H.-R., Yu, S., Oh, J., Kim, K.-H., Lee, J.-H., Moniruzzaman, M., et al. (2019). Nitrate Contamination and Subsequent Hydrogeochemical Processes of Shallow Groundwater in Agro-Livestock Farming Districts in South Korea. Agric. Ecosyst. Environ. 273, 50–61. doi:10.1016/j.agee.2018.12.010

CrossRef Full Text | Google Scholar

Kim, S. H., Kim, H.-R., Yu, S., Kang, H.-J., Hyun, I.-H., Song, Y.-C., et al. (2021). Shift of Nitrate Sources in Groundwater Due to Intensive Livestock Farming on Jeju Island, South Korea: With Emphasis on Legacy Effects on Water Management. Water Res. 191, 116814. doi:10.1016/j.watres.2021.116814

PubMed Abstract | CrossRef Full Text | Google Scholar

Klaus, J., and McDonnell, J. J. (2013). Hydrograph Separation Using Stable Isotopes: Review and Evaluation. J. Hydrol. 505, 47–64. doi:10.1016/j.jhydrol.2013.09.006

CrossRef Full Text | Google Scholar

Lacelle, D., Juneau, V., Pellerin, A., Lauriol, B., and Clark, I. D. (2008). Weathering Regime and Geochemical Conditions in a Polar Desert Environment, Haughton Impact Structure Region, Devon Island, Canada. Can. J. Earth Sci. 45, 1139–1157. doi:10.1139/e08-063

CrossRef Full Text | Google Scholar

Langmuir, D. (1971). The Geochemistry of Some Carbonate Ground Waters in Central Pennsylvania. Geochimica et Cosmochimica Acta. 35, 1023–1045. doi:10.1016/0016-7037(71)90019-6

CrossRef Full Text | Google Scholar

Lee, B.-J., Koo, M.-H., Park, Y.-S., Koh, G.-W., and Park, K.-H. (2006). Hydraulic Diffusivity and Possibility of Conduit-Flow of Groundwater in Eastern Part of Jeju Island. J. Geol. Soc. Korea 42 (3), 439–454.

Google Scholar

Lee, E. S., and Krothe, N. C. (2001). A Four-Component Mixing Model for Water in a Karst Terrain in South-Central Indiana, USA. Using Solute Concentration and Stable Isotopes as Tracers. Chem. Geology. 179, 129–143. doi:10.1016/s0009-2541(01)00319-9

CrossRef Full Text | Google Scholar

Liu, Z., Li, Q., Sun, H., and Wang, J. (2007). Seasonal, Diurnal and Storm-Scale Hydrochemical Variations of Typical Epikarst Springs in Subtropical Karst Areas of SW China: Soil CO2 and Dilution Effects. J. Hydrol. 337, 207–223. doi:10.1016/j.jhydrol.2007.01.034

CrossRef Full Text | Google Scholar

Lorette, G., Lastennet, R., Peyraube, N., and Denis, A. (2018). Groundwater-Flow Characterization in a Multilayered Karst Aquifer on the Edge of a Sedimentary basin in Western France. J. Hydrol. 566, 137–149. doi:10.1016/j.jhydrol.2018.09.017

CrossRef Full Text | Google Scholar

Marfia, A. M., Krishnamurthy, R. V., Atekwana, E. A., and Panton, W. F. (2004). Isotopic and Geochemical Evolution of Ground and Surface Waters in a Karst Dominated Geological Setting: a Case Study from Belize, Central America. Appl. Geochem. 19, 937–946. doi:10.1016/j.apgeochem.2003.10.013

CrossRef Full Text | Google Scholar

Milanović, P. (2021). Dams and Reservoirs in Karst? Keep Away or Accept the Challenges. Hydrogeol. J. 29, 89–100.

Google Scholar

Mocior, E., Rzonca, B., Siwek, J., Plenzler, J., Placzkowska, E., Dabek, N., et al. (2015). Determinants of the Distribution of Springs in the Upper Part of a Flysch Ridge in the Bieszczady Mountains in Southeastern Poland. Episodes. 38 (1), 21–30. doi:10.18814/epiiugs/2015/v38i1/003

CrossRef Full Text | Google Scholar

Moore, P. J., Martin, J. B., and Screaton, E. J. (2009). Geochemical and Statistical Evidence of Recharge, Mixing, and Controls on Spring Discharge in an Eogenetic Karst Aquifer. J. Hydrol. 376, 443–455. doi:10.1016/j.jhydrol.2009.07.052

CrossRef Full Text | Google Scholar

Mudarra, M., Andreo, B., and Mudry, J. (2012). Monitoring Groundwater in the Discharge Area of a Complex Karst Aquifer to Assess the Role of the Saturated and Unsaturated Zones. Environ. Earth Sci. 65, 2321–2336. doi:10.1007/s12665-011-1032-x

CrossRef Full Text | Google Scholar

Palcsu, L., Gessert, A., Túri, M., Kovács, A., Futó, I., Orsovszki, J., et al. (2021). Long-Term Time Series of Environmental Tracers Reveal Recharge and Discharge Conditions in Shallow Karst Aquifers in Hungary and Slovakia. J. Hydrol. Reg. Stud. 36, 100858. doi:10.1016/j.ejrh.2021.100858

CrossRef Full Text | Google Scholar

Parise, M., Ravbar, N., Živanović, V., Mikszewski, A., Kresic, N., Mádl-Szőnyi, J., et al. (2015). “Hazards in Karst and Managing Water Resources Quality,” in Karst Aquifers – Characterization and Engineering. Editor Z. Stevanovic (Springer Series Professional Practice in Earth Sciences), 601–687. doi:10.1007/978-3-319-12850-4_17

CrossRef Full Text | Google Scholar

Park, Y.-Y., Lee, J.-Y., Lim, H.-G., and Park, Y.-C. (2011). Characteristics of the Hydraulic Conductivity of Carbonate Aquifers in Gangwon Province. J. Eng. Geology. 21, 79–85. doi:10.9720/kseg.2011.21.1.079

CrossRef Full Text | Google Scholar

Parkhurst, D. L., and Appelo, C. A. J. (2013). “Description of Input and Examples for PHREEQC Version 3: a Computer Program for Speciation, Batch-Reaction, One-Dimensional Transport, and Inverse Geochemical Calculations,” in Chapter 43 of Section A: Groundwater in Book 6 Modeling Techniques. doi:10.3133/tm6a43

CrossRef Full Text | Google Scholar

Perrin, J., Jeannin, P.-Y., and Zwahlen, F. (2003). Epikarst Storage in a Karst Aquifer: a Conceptual Model Based on Isotopic Data, Milandre Test Site, Switzerland. J. Hydrol. 279, 106–124. doi:10.1016/s0022-1694(03)00171-9

CrossRef Full Text | Google Scholar

Peyraube, N., Lastennet, R., and Denis, A. (2012). Geochemical Evolution of Groundwater in the Unsaturated Zone of a Karstic Massif, Using the Relationship. J. Hydrol. 430-431, 13–24. doi:10.1016/j.jhydrol.2012.01.033

CrossRef Full Text | Google Scholar

Pogačnik, Z., Ogata, K., Pini, G. A., and Tunis, G. (2017). “Karst Structures in Heterogeneous Lithological Units as a Potential Geo-Engineering hazard Factor for Mining and Civil Infrastructures,” in Proceedings of the 3rd Regional Symposium on Landslides in the Adriatic Balkan Region, Ljubljana, Slovenia, October 2017, 11–13.

Google Scholar

Ree, J.-H., Han, R., and Kim, J. K. (2009). Prelithified Deformation of the Jurassic Bansong Group along the Gongsuwon Thrust, South Korea. Geosci. J. 13, 167–173. doi:10.1007/s12303-009-0016-2

CrossRef Full Text | Google Scholar

Rusjan, S., Sapač, K., Petrič, M., Lojen, S., and Bezak, N. (2019). Identifying the Hydrological Behavior of a Complex Karst System Using Stable Isotopes. J. Hydrol. 577, 123956. doi:10.1016/j.jhydrol.2019.123956

CrossRef Full Text | Google Scholar

Schilling, K. E., and Helmers, M. (2008). Tile Drainage as Karst: Conduit Flow and Diffuse Flow in a Tile-Drained Watershed. J. Hydrol. 349, 291–301. doi:10.1016/j.jhydrol.2007.11.014

CrossRef Full Text | Google Scholar

Seo, K. S. (1997). Conodont Fauna from the Gousung limestone, Tanyang Area, Korea and its Biostratigraphic Significance. J. Geol. Soc. Korea. 33, 220–233. (in Korean with English abstract).

Google Scholar

Serno, S., Flude, S., Johnson, G., Mayer, B., Karolytė, R., Haszeldine, R. S., et al. (2017). Oxygen Isotopes as a Tool to Quantify Reservoir-Scale CO 2 Pore-Space Saturation. Int. J. Greenhouse Gas Control. 63, 370–385. doi:10.1016/j.ijggc.2017.06.009

CrossRef Full Text | Google Scholar

Somers, L. D., and McKenzie, J. M. (2020). A Review of Groundwater in High Mountain Environments. WIREs Water. 7, e1475. doi:10.1002/wat2.1475

CrossRef Full Text | Google Scholar

Stueber, A. M., and Criss, R. E. (2005). Origin and Transport of Dissolved Chemicals in a Karst Watershed, Southwestern Illinois. J. Am. Water Resour. Assoc. 41, 267–290. doi:10.1111/j.1752-1688.2005.tb03734.x

CrossRef Full Text | Google Scholar

Tobin, B., and Schwartz, B. (2012). Quantifying Concentrated and Diffuse Recharge in Two marble Karst Aquifers: Big Spring and Tufa Spring, Sequoia and Kings Canyon National Parks, California, USA. Jcks. 74 (2), 186–196. doi:10.4311/2011jcks0210

CrossRef Full Text | Google Scholar

Tritz, S., Guinot, V., and Jourde, H. (2011). Modelling the Behaviour of a Karst System Catchment Using Non-Linear Hysteretic Conceptual Model. J. Hydrol. 397, 250–262. doi:10.1016/j.jhydrol.2010.12.001

CrossRef Full Text | Google Scholar

Vesper, D. J., and White, W. B. (2004). Storm Pulse Chemographs of Saturation Index and Carbon Dioxide Pressure Implications for Shifting Recharge Sources During Storm Events in the Karst Aquifer at Fort Campbell, Kentucky/Tennessee, USA. Hydrogeol. J. 12, 135–143. doi:10.1007/s10040-003-0299-8

CrossRef Full Text | Google Scholar

Waltham, T., and Lu, Z. (2007). “Natural and Anthropogenic Rock Collapse over Open Caves, Geological Society, London, Special Publications,” in Natural and Anthropogenic Hazards in Karst Areas: Recognition, Analysis and Mitigation. Editors M. Parise, and J. Gunn (Geological Society London, Special Publications), 279, 13–21. doi:10.1144/SP279.3

CrossRef Full Text | Google Scholar

Zhu, H., Xing, L., Meng, Q., Xing, X., Peng, Y., Li, C., et al. (2020). Water Recharge of Jinan Karst Springs, Shandong, China. Water. 12, 694. doi:10.3390/w12030694

CrossRef Full Text | Google Scholar