Hydrogeochemical Characteristics and Processes of Shallow Groundwater in the Yellow River Delta, China

The Yellow River Delta is one of the biggest river deltas in China, and the shallow groundwater plays an important role in the development of the local agriculture and ecosystem. However, people are still unclear about the hydrochemical characteristics and mechanisms of the shallow groundwater. In this study, the authors collected and analyzed 81 groundwater samples from the delta plain and piedmont alluvial plain, and explored the hydrochemical features and causes through Piper diagrams, correlation analysis, ionic ratios, and speciation calculations. The results showed that anions were dominated by Cl and HCO3, the concentration of which was much more than that of SO4 and CO3. The groundwater can be divided into various types, including Na–Cl, Ca–Mg–HCO3, Na–HCO3 and Ca–Mg–Cl. This study tested an alternative method–ionic ratios based on the cumulative frequency distribution to characterizing the hydrochemical groups. According to different ion ratios and hydrogeological conditions, three hydrogeochemical zones with different dominant factors have been determined: Weathering—Fresh Water Zone (Zone I), Evaporation—Saline Water Zone (Zone II), and Seawater Mixing Zone (Zone III). As the calculated saturation index show, the calcite and dolomite are saturated, while the halite and gypsum from Zone I to Zone III tend to be saturated. In addition, cation exchange is an important hydrochemical process in the area, and Zone III experiences inverse ironic exchange. In conclusion, this hydrogeochemical zonation would be favorable for water resource management in the Yellow River Delta.


Introduction
Groundwater is an important resource for domestic, agricultural, industrial, and ecological development in arid and semi-arid regions [1]. Exploring groundwater chemistry characteristics, origin, and evolution are conductive to understanding groundwater circulation, hydrodynamic conditions and environmental change, and can better regulate the management and protection of local groundwater resources [2,3]. Extensive research has been conducted on the hydrochemical evolution and suitability of groundwater in basin and river delta [4][5][6]. These studies suggested that hydrogeochemical processes of shallow groundwater are complicated and are influenced by climate, sedimentary environment, and anthropogenic activities. The hydrogeochemical processes of shallow groundwater include dissolution or precipitation, evaporation, ion exchange, and mixing [7][8][9]. The mathematical statistics, multivariate analysis, ionic ratio, mineral saturation index, inverse modeling, and isotopic tracing are regarded as the major methods for studying groundwater chemistry [10][11][12][13].
Located in the southwest coast of Bohai Sea, the Yellow River Delta (YRD) is one of the largest river deltas in China. With rich oil and gas and land resources, the Yellow River Delta is an important energy base and agricultural area [14]. The exploitation of the shallow underground brackish water and salt water has become a major way to deal with the shortage of agricultural water resources in this area [15]. Featured by large loose-sediment content and frequent siltation, burst, and migration, the Lower Paleo-Yellow River caused complex stratigraphic texture [16]. The ecosystems of the Yellow River Delta severely degrade, which attributes to the salinization caused mainly by seawater intrusion [17]. In addition, the buried depth of the shallow groundwater table in the region is very shallow and may experience strong evaporation of groundwater. However, the previous studies mainly focused on the influence of seawater incursion on the shallow groundwater environment [17][18][19], and there are a limited number of studies on the hydrochemical features and the formation mechanism of groundwater [20].
Therefore, this study aims to: (1) investigate the hydrochemical characteristics of the shallow groundwater; (2) identify different hydrochemical zones, revealing the predominant hydrogeochemical processes. The research results aim to provide helpful information for the utilization and protection of shallow groundwater resources in the Yellow River Delta.

Study Area
The Yellow River Delta is located at the mouth of the Yellow River and borders the Bohai Sea to the east and north ( Figure 1). Topographically, the YRD can be divided into the alluvial-marine plain in the coastal areas and the piedmont proluvial-alluvial plain situated in the north of Shandong hills [21]. The middle reaches of the Yellow River carry a large amount of loose sediment to downstream through the Loess Plateau. The Yellow River transports 1.08 billion tons of sediment annually, and one-fifth of the sediment is silted near the entrance [22]. On the south of the study area, coarse sediments come from rock wearing of the Shandong hills located to the south of the study area. The YRD has a warm temperate semi-humid monsoon climate with an average annual rainfall of 530-630 mm, which is far less than the average evaporation of 1470-2246 mm.
Originating from the terrigenous and marine sources in the delta plain, the Quaternary sediments have a thickness of over 400 meters. The western coast of Bohai Bay went through multiple Quaternary marine transgressions induced by sea-level changes during the glacial and interglacial eras [23]. The Quaternary aquifer is divided into three groups, including unconfined and semi-confined aquifers (Aquifer I), middle confined porous aquifers (Aquifer II), and deep confined porous-fissure aquifers (Aquifer III) [24]. Aquifer I is mainly composed of Holocene and late Pleistocene fine sand, silt and clay silt. The buried depth is less than 60 m. The lower Pleistocene semi-confined aquifers are hydraulically connected with the Holocene unconfined aquifers, so the Aquifer I can be effectively categorized into a single shallow aquifer. The depth of the shallow groundwater table in the delta plain ranges from 0.5 to 2.5 m under the surface. Precipitation and river infiltration are the main groundwater recharge ways, and evaporation is the main way of groundwater discharge. The groundwater flows from southwest to northeast with a low velocity.
In the southern part of the study area, there exposes a large area of Cambrian-Ordovician limestone and Archean granite. Originating from bedrock wearing, the aquifer medium in the alluvial-proluvial plain is composed of fine sand and medium-fine sand, with medium-coarse sand in some regions. The depth of the shallow groundwater table is generally between 1.8 and 17 m under the land surface. This area is rich in shallow freshwater resources, and artificial exploitation can be mainly used for discharging groundwater. Moreover, the long-term over-exploitation of groundwater has continuously decreased the groundwater level and caused saltwater intrusion. The groundwater mainly flows from the southwest to the northeast.

Groundwater Sampling and Laboratory Analysis
In the study area, eighty-one groundwater samples were taken from irrigation wells with a depth of 20 to 60 m, in September 2015 ( Figure 1). The well was purged at least half an hour earlier before sampling so as to assure that the groundwater samples were fresh. A WTW-Multi-340i/SET (Germany) multi-parameter measuring instrument was used to measure the groundwater temperature, pH value and redox potential (Eh) on site. Groundwater samples were collected at 0.5 m below the water table that then went through 0.45 µm membrane filters. Groundwater samples for major cations analyses were acidified with HNO 3 to pH < 2, and the unacidified samples were used for the major anion analysis. All samples were stored in portable ice-filled coolers at 4 • C during transport and preserved in refrigerator at the temperature of 4 • C in the laboratory until analysis.
All the chemical component tests were carried out at the Water Quality Analysis Laboratory of the Institute of Hydrogeology and Environmental Geology. Major cations (Na + , K + , Ca 2+ , and Mg 2+ ) were determined through inductively coupled plasma-mass spectrometry (Agilent 7500ce ICP-MS). Ion chromatography (DX-120, Dionex) was adopted for determining major anions (Cl − , SO 4 2− , and NO 3 − ). Concentrations of HCO 3 − were measured through acid-base titration. Total dissolved solids (TDS) were determined by gravimetric analysis. The analytical precision error of the ICP-MS and DX-120 instruments were both less than 5%. Among most groundwater samples, ion charge imbalances were <5%. Descriptive statistics of chemical components are shown in Table S1.

Data Analysis
The groundwater chemistry diagram (Piper) was analyzed by RockWare AQQA (1.1) software. The ARCGIS (10.0) was applied to map the spatial location of water chemical ions. The Pearson correlation coefficient used to assess the correlation coefficient between hydrochemical indicators was calculated by IBM SPSS Statistics (19.0). The OriginPro 2016 was used for general data analysis and graphing, such as the cumulative frequency distribution curve, Gibbs diagram, and ionic exchange. In addition, mineral saturation index was calculated by geochemical simulation program PHREEQC (Version 3.3.5) [25].

General Groundwater Chemistry
The pH values of groundwater ranged from 6.8 to 8.2 (average 7.5), so it indicates that the environment is of the neutral pH value (Table S1). The TDS concentrations varied widely from 328 mg/L to 145,997 mg/L, with an average of 15,027 mg/L ( Figure 2). The TDS concentration in delta plain was abnormally high, while that in piedmont plain was relatively low ( Figure 3). Typically, in groundwater samples throughout the whole study area, the concentration is Na > Mg > Ca > K (n = 39) or Na > Ca > Mg > K (n = 21) ( Figure 2). Other groups of cations were found only in a few sites across the delta. Anions are dominated by Cl and HCO 3 , the concentration of which was much more than that of SO 4 and CO 3 . The Piper diagram ( Figure 4) shows different groundwater types, including Na-Cl, Ca-Mg-HCO 3 , Na-HCO 3 and Ca-Mg-Cl types. The Na-Cl type dominates the groundwater in the delta plain, while the groundwater in the alluvial-proluvial plain in the southern part of the study area is mainly dominated by Ca-Mg-HCO 3 type. A few groundwater samples were of Na-HCO 3 and Ca-Mg-Cl types.  Table 1 shows the Pearson correlation coefficient of the hydrochemical parameters. There was a positive correlation among TDS, Na + , Mg 2+ , SO 4 2-, and Cl -, with the correlation coefficients greater than 0.9 (p < 0.01). It indicates that the groundwater is mainly composed of Na + , Mg 2+ , SO 4 2-, and Clin terms of the hydrochemical composition, and these components may have the same sources. The relationships between HCO 3 and other chemical parameters are weak (coefficient < 0.62, p < 0.01), which suggests that HCO 3 may originate from different hydrogeochemical processes, compared with other chemical parameters.

Ionic Ratios
Ionic ratios can be effectively used to reveal hydrochemical characteristics of groundwater [9,26]. In the present study, Ca/Cl, Mg/Ca and Na/HCO 3 ionic ratios were assessed based on cumulative frequency distribution diagram so as to accomplish hydrochemical groupings. The results of these tests are shown in Figures 5-7. The cumulative frequency distribution curve for Ca/Cl ratio (Figure 5a) can be used to divide the samples into four groups and the distribution of these groups on the Piper diagram are distinctly separated. Group A1 (0-0.1) samples, accounting for 28.40% of total samples, take relatively higher chloride proportions. Therefore, they are mainly located on the far right of the Piper diagram (Figure 5b), surrounding seawater concentration. Based on Figure 5c, samples of Group A1 were mainly distributed in the coastal areas, and the water chemical composition is obviously affected by seawater mixing. Group A2 (0.1-0.4) samples occupied 25.9% of the dataset. Although there were two samples of the Na-HCO 3 type, most samples were of the Cl-Na type waters situated on the left of Group A1 samples in the Piper diagram (Figure 5b). Geographically, these samples were distributed on the west of group A1 samples, so it means that the seawater mixing is relatively weak. Groups A3 (0.4-1.4) and Group A4 (>1.4) represented 33.3% and 12.3% of the whole samples respectively. Groups A3 samples mainly took the central part of the diamond field where HCO 3 was dominant over Cl -, and a few samples were placed in the centre where SO 4 2was the dominant anion. The dominant cation in Group A3 did not take a prominent role, with Mg 2+ and Ca 2+ ranking the first and second places, followed by Na + . Group A4 was clearly clustered in the left-hand side of the Piper plot. These samples could illustrate the conditions in the proluvial-alluvial plain, and could represent higher Ca 2+ and HCO 3 concentrations produced by strong leaching. Ca and Mg are two common divalent cations in natural groundwater. As seawater has higher content of Mg 2+ content than the Ca 2+ , Mg/Ca ratio is often used to evaluate the degree of seawater intrusion [27,28]. The cumulative frequency distribution curve for Mg/Ca ratio (Figure 6a) shows that there are four groups as well. Group B1 samples (0-1.0) occupy about 19.8% of the total samples and form two independent clusters. The first cluster is dominated by Ca-HCO 3 with low Mg and Na content (Figure 6b). The second cluster is a tight grouping of Ca-Cl type, which may be related to cation inverse changes [29]. The samples are mainly distributed in the piedmont plain, and a few are distributed along the Yellow River. The groundwater in this area is mainly recharged by precipitation and surface water infiltration. Group B2 (1.0-2.0) samples are evenly distributed across the diamond-field on the Piper diagram (Figure 6b). The wide distribution of Group B2 samples in the study area (Figure 6c) suggests this group is not defined by a specific area and probably reflect a background population where Na-Cl and Ca-Mg HCO 3 type groundwater are most common. This group accounts for 42.0% of the total dataset. Group B3 (2.0-4.0) and Group B4 (>4.0) which comprises 38.3% of samples, which is obviously concentrated on the right-hand side of the diamond-field on the Piper diagram. Geographically, these samples are mostly confined to the coastal zone of the alluvial-marine aquifer, showing high Mg concentrations by seawater mixing. Similarly, Na/HCO 3 ratio reflects a contrast of solutes typically derived from terrestrial and marine sources. Group C1, C2, C3 and C4 occupy 33.3%, 30.9%, 25.9% and 9.9% of the total dataset. Seen from Piper diagram (Figure 7b), these groups are sequentially distributed across the diamond-field from left to right as the scale of Na to HCO 3 variation. The distribution patterns of these groups were similar to other groups already described (Figure 7c). Group C1 represents a zone where groundwater is recharged by fresh HCO 3 water. Group C2 is particularly distributed in the transition zone between piedmont plain and delta plain, which reflects the background group of the sample and has no obvious ionic effect. Instead, this group may be more indicative of combined processes, including dissolution, ionic exchanges, or mixing of other groups. In addition, the higher Na contents (Group C3 and C4) mainly occurs in the alluvial-marine aquifers in which seawater has significant influences on groundwater composition.

Hydrogeochemical Zonation
The result shows that evaluating the probability distribution of ion relations can help distinguish variations in ionic relationships, and different hydrochemical groups can be defined by critical points. For the purpose of revealing the complicated hydrogeochemical processes for controlling groundwater chemistry, Ca/Cl, Na/Cl and Mg/Ca ratios and hydrogeological settings are used to identify three hydrogeochemical zones (Figure 8).

Weathering-Fresh Water Zone (Zone I):
This zone is mainly situated in the piedmont proluvial-alluvial plain and is the groundwater recharge zone. The aquifer sediments are mainly composed of fine sand and medium sand with relatively high permeability, and the hydraulic gradient of about 5‰ [30]. The strong dissolution brings about low TDS concentrations (average of 100.5 mg/L), high concentrations of HCO 3 -(average 429.3 mg/L), and relatively low concentrations of Cl -(219.7 mg/L) (Table S2). It is characterized by high Ca/Cl (>0.4) and low Mg/Ca (<1.0) and Na/HCO3 (<1.0). On the Gibbs diagram, groundwater samples in Zone I are mainly located in the rock weathering dominance (Figure 9a), so it indicates that the hydrochemical composition are mainly controlled by water-rock interaction. The TDS concentration of the three samples in this zone is abnormally high (2050-4200 mg/L), which was caused by the recharge of salt water in the east and north after the excessive exploitation of groundwater. [31].
Evaporation-Saline Water Zone (Zone II): From Zone I to Zone II, the aquifer medium varies from coarse to fine particles, and the permeability is relatively low [30]. According to the survey, the hydraulic gradient of the shallow groundwater in this zone is <0.3‰, and buried depth of water table ranges between 0.5-2.5 m. As a result, groundwater evaporation and concentration in this area is strong. As it can be seen from Figure 9a, the groundwater samples in Zone II lie between the rock weathering dominance and seawater instruction dominance, which means that the evaporation and concentration are important hydrochemical processes of groundwater composition. The TDS concentration in Zone II is significantly higher than that in Zone I. The average TDS concentration in this zone is 3839 mg/L. Na + is the main cation with an average of 1003.7 mg/L, followed by Ca 2+ (average 167.0 mg/L) and Mg 2+ (average 168.4 mg/L) (Figure 8b). The dominant anion is Cl -, followed by HCO 3 and SO 4 2-, with average values of 1669.7 mg/L, 664.9 mg/L and 428.9 mg/L, respectively (Figure 8b). Seawater Mixing Zone (Zone III): This zone is mainly located in the coastal areas. Since the Late Pleistocene, coastal areas have experienced three major transgressions in the western Bohai Sea [32], so the massive paleo-seawater occurs in the alluvial-marine aquifers. In addition, modern seawater intrusion is usually caused by upward tracking of seawater, storm surges, and over-exploitation of groundwater [19]. Groundwaters sampled in Zone III mainly come from the seawater dominance (Figure 9a), which shows that seawater mixing will have a significant impact on the chemical composition of groundwater. The chemical composition of groundwater is characterized by lowest Ca/Cl (<0.1) and highest Mg/Ca (>4.0) and Na/HCO3 (>6.3). After long-term evaporation and concentration, the saline water in the marine sediments could become brine with exceedingly high TDS [33]. Underground brines are widely distributed in coastal areas in the western Bohai Bay, and its buried depth is less than 60m [34]. In Zone III, the maximum concentration of TDS reaches 146.0 g/L with an average of 40.6 g/L. Compared with Zone I and Zone II, the concentration of Mg 2+ is significantly higher than Ca 2+ , and that of SO 4 2is significantly higher than HCO 3 in groundwater of Zone III (Figure 8b), which is notably similar to the chemical component of seawater [35].

Hydrogeochemical Processes
As discussed above, the shallow groundwater chemistry is mainly affected by rock weathering, evaporation, and seawater mixing. In order to determine the sources of solutes in different zones of groundwater, this study drew several bivariate plots. The scatter plots of (Mg 2+ /Na + ) versus (Ca 2+ /Na + ) and (HCO 3 − +Na + ) versus (Ca 2+ +Na + ) were used to assess the relative contribution of the three main weathering processes, including carbonate, silicate, and evaporite weathering [36,37]. Based on Figure 9b-c, Zone I groundwater samples are located between carbonate and evaporate weathering dominance, which indicates that groundwater chemical components in Zone I mainly come from dissolution of carbonate minerals and silicate minerals in southern bedrock area. The bivariate plots show most of the samples in Zone II are correlated with silicate weathering dominance (Figure 9b-c). In contrast, the Zone III samples approach the evaporate dominance. The significant correlation (r = 0.987) between Na + and Clindicates dissolution of evaporate minerals, such as halite.
Mineral equilibrium calculations are conductive to determine the prevailing thermodynamic processes in groundwater system [38]. Thus, saturation indices were used to further determine the saturation and precipitation (over-saturation) and dissolution (under-saturation) of the possible mineral phases. SI values were calculated using geochemical simulation program PHREEQC (Interactive 3.3.5) and phreeqc database. As shown in Figure 10a,b, SI values for calcite and dolomite of all the groundwater samples were greater than zero (SI > 0), and SI values of Zone II and Zone III groundwater samples were greater than that of Zone I. This phenomenon indicates that groundwater is super-saturated with calcite and dolomite, and these minerals are prone to precipitate. In addition, SI values for gypsum and halite of all the groundwater samples were less than zero (SI > 0) (Figure 10c,d). From Zone I to Zone III, SI values of gypsum and halite gradually increased as salt content increases, which indicates that these minerals trend to be saturated. In order to assess the influence of cation exchange and adsorption reactions on groundwater chemistry of the study area, a scatter plot of (Ca 2+ +Mg 2+ -HCO 3 --SO 4 2-) versus (Na + +K + -Cl -) was made (Figure 9d). All the groundwater samples were close to the y = − x line, which indicates that the cation exchange probably has a significant influence on groundwater chemical composition of the Yellow River Delta. In Zone III, groundwater samples were located in the second quadrant, which shows a reverse cation exchange with increasing trend of Ca 2+ and decreasing trend of Na + . In many coastal areas, such as Western Senegal [36], southeastern South Korea [26] and southern Greece [37], there have been reports about reverse cation exchange associated with seawater intrusion. In Zone I and Zone II, the cation exchange capacities were lower and the reaction direction was complicated, compared with Zone III (Figure 9d).

Conclusions
In the study, TDS concentrations varied widely, ranging from 328 to 145,997 mg/L. Concentrations of major anions were Cl > HCO 3 > SO 4 > CO 3 , while the major cations were Na > Mg > Ca > K. Based on Ca/Cl, Mg/Ca and Na/HCO 3 ionic ratios, Piper diagrams and hydrogeological conditions of the region, three distinctive hydrogeochemical zones were determined. Zone I is located in the floodplain in the piedmont of the recharge area. The TDS concentration of groundwater is low, and the groundwater is mainly Ca-Mg-HCO 3 type. Rock wearing is the predominant process for the formation of groundwater components. Zone II is situated in the intermediate between the proluvial-alluvial plain and the delta plain. Under stagnant water conditions, the hydrochemistry is mainly affected by evaporation and concentration. With the increase of salinity, the groundwater correspondingly evolves into Na-Cl type. Zone III is located in the delta plain around the coast. The local hydrochemistry is closely related to seawater mixing and is characterized by abnormally high TDS. In addition, cation exchange reaction was also a significant hydrogeochemical process in the study area. This study provides a feasible means for discriminating hydrochemical characteristics and estimating hydrogeochemical processes. Due to the lack of groundwater level dynamic data, seawater intrusion simulation was not performed in this study. In future research, we would monitor groundwater level and develop field tests for key hydrogeological parameters to support hydraulic simulation.
Supplementary Materials: The following are available online at https://www.mdpi.com/2073-444 1/13/4/534/s1. Table S1, Statistics of major chemical parameters of groundwater in the study area (n = 81); Table S2, Statistics of major chemical parameters of groundwater in different zones.