Abstract
A complete understanding of groundwater circulation as well as the transport and distribution of solutes in arid-semiarid basin regions is a prerequisite for the safe use of groundwater resources. The distributions of the groundwater flow systema and solutes are affected by the basin morphology, lithology, and climate variations; therefore, they can change over geologic time. In this study, we performed a case study of the Qaidam Basin in the northeastern Tibetan Plateau, in which we utilized reactive solute transport simulations to build a numerical model in TOUGHREACT for a typical section of the eastern Qaidam Basin since the last glacial period. The results show that the groundwater in the eastern Qaidam Basin developed into a three-level groundwater flow system and that the seepage velocity of the local water flow system is significantly higher than that of the intermediate and regional water flow systems. Although groundwater in the discharge region has been continuously concentrated and enriched since the last glacial period, the distributions of the groundwater flow system and solutes have been greatly affected by climate variations. During warm periods, the centres of groundwater discharge and solute concentration shifted to areas with more groundwater recharge; in contrast, both centres shifted to the central basin during drought periods. The groundwater in the basin mainly contains Na+ and Cl− ions, which vary significantly from the recharge region to the discharge region. Evaporation of groundwater results in increases in the concentrations of most of the components except HCO3−. The groundwater in the discharge region is currently in the stage of carbonate precipitation and is far from gypsum and halite precipitation.
1. Introduction
As an important component of water resources, groundwater is the main and sometimes only source of water in the arid and semiarid regions of China. More than 50% of the water supply of the major cities in northern and western China comes from the groundwater, and this percentage is as high as 80% in many cities [1]. The reasonable development and utilization of groundwater resources has become the main factor affecting the sustainable development of local economies, and water security is a significant issue that China urgently needs to address [2].
The Qaidam Basin is located in the arid and semiarid regions of China. It is one of the four major inland basins of China and contains abundant mineral resources, such as salt and borax. However, water resources are extremely scarce, which will seriously threaten the sustainable development of the economy and the ecological security of the Qaidam Basin. The development and utilization of water resources has not been prudent, which has caused environmental and geological problems, such as an expanding groundwater depression cone, saline water, and vegetation degradation [3]. In addition, in recent years, global climate warming has significantly influenced the Qaidam Basin, making it the region experiencing the most significant climate changes in the Tibetan Plateau [4]. These changes have also considerably influenced the local water resources [4–6]. Therefore, a comprehensive understanding of the groundwater flow in the Qaidam Basin is key to the reasonable development and utilization of local water resources.
The movement and evolution of groundwater in a basin is controlled and affected by multiple factors. Tóth summarized two major classes of controlling factors: the geometric morphology of the basin and the geological conditions of the basin [7]. In addition, many scholars believe that the upper boundary conditions and the physical properties of the fluid are also important factors that affect the movement of groundwater in a basin [8–11]. Li and Hao [12] adopted a multidisciplinary research method to determine that inland basins can be roughly divided into four grades of flow systems and three characteristic zones of salt migration, and they explained the formation mechanism for the freshwater segment of an inland basin. Lin and Jin [13] began with the principles of the hydrologic cycle and considered the motion of water and salt as one entity to systematically study the characteristics of hydrologic circulation and salt migration and accumulation on different levels under the natural and human-influenced conditions of a basin.
Many studies have focused on the groundwater flow and saltwater migration in the Qaidam Basin [14–21]. For example, Ma et al. [19] performed an isotope study in the Golmud subdrainage basin and concluded that groundwater was the main source of the salt lakes. Ye et al. [21] investigated the hydrochemical characteristics and sources of brine in the Gasikule salt lake in the northwestern Qaidam Basin and concluded that the sources of the salt in the lake included stream water, the leaching of Pliocene salt-bearing host rocks, and Ca-Cl-type deep water. Chen et al. [14] used rare earth elements as tracers to study the dissolution of source rocks in different water bodies in the Qaidam Basin. Hu and Jiao [15] constructed a regional three-dimensional groundwater flow model for the Qaidam Basin that was calibrated using GRACE data from 2003 to 2012. However, these studies mainly focused on short time scales. Few studies have explored the groundwater flow and salt transport over long time scales, such as the 1000-year or 10,000-year scales. In addition, previous studies mainly focused on the western part of the Qaidam Basin, where the hydrogeologic conditions are different from those in the eastern part of the basin [22].
In this paper, we use TOUGHREACT [23] to perform a reactive solute transport simulation and consider the influence of water-soluble components on the fluid density. We take a typical cross section of the Qaidam Basin as an example and simulate the migration and evolution of the groundwater and salt in the eastern Qaidam Basin since the last glacial period. Based on the simulation results, we identify multiple levels of the groundwater flow system and present the distribution of the solutes in the groundwater during different geological periods. This information is important for the sustainable utilization of local groundwater resources and for coping with climate change.
2. Study Area
2.1. Geological and Hydrogeological Setting
The Qaidam Basin is a large closed intermountain fault basin in which the maximum thickness of the Quaternary sediments exceeds 3000 m [11]. From the piedmont to the central basin, the structural characteristics of the water-bearing media generally have zonal or semizonal distributions. In the vertical direction, from the piedmont to the central basin, the water-bearing unit changes from a single phreatic aquifer to a multilayer confined aquifer. The lithology of the aquifers changes from coarse to fine-grained, and the thickness decreases. The water yield and permeability change from high to low.
Piedmont alluvial plains are distributed around the basin. The aquifers in these plains are dominated by sand and gravel [22]. The aquifer thickness varies from tens of metres to one to two hundred metres, and the hydraulic conductivity varies from several metres per day to hundreds of metres per day [22]. These areas are the water-rich regions of the basin, and the daily single-well water yields vary from several tens of cubic metres to three thousand cubic metres [22]. Alluvial plains are located in the central basin. The aquifer changes from a single-layer structure to a multilayer structure, and the groundwater transitions from phreatic water to confined water. The aquifer is generally no more than 50 m thick, and the thinnest part is less than 10 m thick [22]. The lithology of the aquifer is generally moderately fine sand, coarse sand, and fine sand, and the aquiclude layers are mainly composed of impermeable sandy loam, loam, and clay. The hydraulic conductivity of the aquifer varies relatively significantly, which is generally less than one metre per day; however, the maximum value is tens of metres per day [22].
The climate of the Qaidam Basin is classified as an extremely arid inland climate. The basin is surrounded by high mountains (Figure 1), and the flow of warm and moist air is restricted. Therefore, precipitation is rare, and potential evaporation is significant. Based on an analysis of multiyear observational data at meteorological stations, the average rainfall in the basin is 16.09–189.73 mm, and the average potential evaporation in the plains area is 1973.62–3183.04 mm [17]. The potential evaporation rates are much higher than the precipitation rates; therefore, rainfall only slightly replenishes the groundwater. The groundwater mainly receives vertical infiltration from rivers in the alluvial fans as well as lateral supply from channel seepage, infiltration of irrigation water, and water from bedrock fissures [22]. The seepage conditions in the alluvial fans are good. The groundwater rapidly moves to the overflow zone at the front edge of the alluvial fans, and some of the groundwater flows to the surface in the groundwater spill belt to form rivers. A portion of the groundwater continuously seeps towards the central basin and eventually evaporates from phreatic water or is directly recharged to the terminal lake.
The distribution of chemical components in the groundwater in the basin is mainly controlled by the groundwater recharge and discharge conditions. The chemical components of the groundwater in the middle and upper parts of the alluvial fan are not significantly different from those of the river water. Based on a chemical analysis of well CK7, which was drilled in the middle of the alluvial fan and sampled at several depths in 2013, the average total dissolved solids (TDS) of the groundwater is 403 mg/L. The TDS of the river water in the Nuomuhong mountain pass is 427 mg/L [24]. When the groundwater flows to the overflow zone, the concentrations of the chemical components in the groundwater rapidly increase by evaporation, and the TDS of the groundwater can reach the concentration of saline water. The TDS of the groundwater can exceed the concentration of brine until it reaches the region adjacent to the terminal lake in the central basin.
2.2. Selected Cross Section and Salinization Period
The typical cross section selected in this study is located in the Nuomuhong area in the eastern Qaidam Basin. The cross section begins at the top of the Nuomuhong alluvial fan in the south and extends to the Amunike Mountain piedmont in the north. Because the Nuomuhong alluvial fan has a regular form and the surface runoff is relatively high, it is suitable for studying the variations in the groundwater flow system caused by surface water [24]. From the top to the edge of the alluvial fan, the direction of the cross section is aligned with the slope gradient as well as the groundwater flow direction. In the central basin, the groundwater and surface water flow to the west. However, the groundwater seepage velocity in the east-west direction (perpendicular to the cross section) is extremely slow and does not have a significant influence on the flow simulation along the section [25].
The Nuomuhong area is located in the eastern Qaidam Basin, and it has not entered the salt precipitation stage until recently. Based on the analytical results of the Dacan-1 well on the eastern coast of Dabuxun Lake, the salinization zone and desalination zone from 3200 m to 45 m depth indicate that the Qarhan lake water at that time was in the transition between saline water and brackish water, and the lake was likely dominated by brackish water [26]. The age of the most recent appearance of a turf or humus soil layer inferred from a pore is 70 ka BP, which is close to the start of the last glacial period at 75 ka BP; this indicates that after 70 ka BP, the Qarhan Salt Lake began to develop. Therefore, the salinization period is determined to have extended from the beginning of the last glacial period at 70 ka to the present.
2.3. Paleohydrologic Reconstruction
Rainfall is sparse in the Qaidam Basin, and the river runoff from the mountainous area represents most of the groundwater recharge resources in the basin. Therefore, the key to determining the basin source and sink terms is to determine the river runoff and the potential evaporation. The river runoff is affected by the amount of rainfall in the catchment area, and the rainfall is affected by regional and global climate change. Therefore, we can approximate the rainfall rates under different climate conditions and further estimate the river runoff based on the rainfall rates. According to a hydrogeological survey in the Qaidam Basin, the variation in river runoff in this area is 60% of that of precipitation [22]. The period of 45–25 ka BP was a relatively warm period during the last glacial period, and the temperature reached that of an interglacial climate [27]. The temperature in the southwestern Tibetan Plateau was 3–4°C higher than at present, and the amount of rainfall was 40%–100% higher. The amount of rainfall in the Qinghai Lake region adjacent to the Qaidam Basin was 79% higher. The Qarhan Salt Lake was also three to seven times larger than its current area. The highest temperatures during the Holocene warm period occurred at 7.2–6.0 ka BP. The temperature in the Qinghai Lake area was 3°C higher than the current temperature, and the precipitation was 70%–80% higher [27]. The Golmud River flood event in 1989 caused the Qarhan Salt Lake area to expand by approximately 3.9 times [28]. Based on hydrological monitoring data, the average temperature in 1989 was approximately 2°C higher than the multiyear average. The rainfall in the Golmud River basin during that year was 3.2 times the multiyear average, and the river runoff was 2.1 times the multiyear average [22].
These data show that there is a positive correlation between changes in precipitation and temperature changes in the Qaidam area. Based on these data, we preliminarily infer the paleohydrological characteristics of the Qaidam Basin since 70 ka by combining the characteristics of global climate change from the Late Pleistocene epoch to the Holocene (Figure 2). The current multiyear average temperature in the Qaidam Basin is approximately 4°C. Based on the global climate characteristics, the temperature during the first interglacial period was approximately 7–8°C. According to the results using the isotope and rare gas method in this research, the Pleniglacial temperature was approximately 1–2°C [25]. The first subsidiary glacial stage was warmer than the Pleniglacial stage, and it was colder than the first subsidiary interglacial stage, which is believed to be the same as today. As described above, the variation in precipitation is positively correlated with temperature. We use the current amount of rainfall as a reference to infer the variations in rainfall at other times. The river runoff is positively correlated with the rainfall, but the magnitude of the fluctuations in runoff is smaller than the magnitude of the rainfall variation. Based on the analysis of multiyear data, the magnitude of the variation in runoff in the Golmud region is approximately 60% of that of the rainfall [22]. The evaporation capacity is negatively correlated with both the amount of precipitation and the temperature. A temperature change of 2°C causes the evaporation to fluctuate by approximately 20%; based on this correlation, we can infer the evaporation at different times [22].
3. Simulation Approach
3.1. Conceptual Model
The horizontal distance of the cross section is 63 km. The vertical range of the simulation is from the ground surface to the Quaternary basement. The simulation thickness gradually increases from the top of the alluvial fan to the central basin, and the maximum simulation thickness is 1800 m. The bottom is defined as an impermeable boundary because the exchange between the Quaternary aquifer and the bedrock is minimal and can be ignored [22]. The upper boundary is defined as a flow boundary or Dirichlet boundary. The piedmont river infiltration is a given flow boundary, which is the only supply source in the simulation region. We determine the amounts and locations of infiltration based on the river infiltration distribution. The spring, river, and groundwater evaporation are defined as Dirichlet boundaries for the unidirectional discharge.
From the edge of the basin towards the centre, the sediment is mainly composed of sandy gravel, sand, silty soil, silty clay, and clay. Because limited data are available from in situ pumping tests along the cross section, the hydrogeological parameters of the different lithologies are based on empirical values [22]. The hydrogeological parameters were adjusted during the model calibration process by comparing the simulation results with observation data. The calibrated parameters are shown in Table 1. The discussion of the simulation results is based on the calibrated parameters.
Before the unconsolidated sediment entered the basin, it had undergone sufficient weathering and hydrolysis, and the soluble salt had been fully filtered. After entering the basin, most of the clastic sediments are composed of relatively insoluble minerals that are not easily weathered. A comparison of the chemical analysis of well CK7 with that of the surface water demonstrated that the groundwater in the basin will not dissolve additional components. Therefore, we can ignore the soluble salt formed by mineral weathering in the basin. It is assumed that all of the soluble salt came from outside the basin by surface water recharge.
The movement of water-soluble components in an aquifer mainly includes convection and hydrodynamic dispersion. The hydrodynamic dispersion coefficient is replaced by the diffusion coefficient in TOUGHREACT [23]. The groundwater flow rate is relatively slow in the central basin, so it is reasonable to make the replacement. The value of the diffusion coefficient is set to 2.95 × 10−6 m2/s [24]. The chemical reactions considered in the model mainly include mineral dissolution/precipitation, convective mixing, and cation exchange. The hydrochemical reactions are based on the local equilibrium model of water-soluble components. The ionic activity coefficient of high ionic strength is calculated based on the Pitzer model of ion interaction [23, 29].
The main ion components in the groundwater are the main variables in the hydrochemical model, including K+, Na+, Ca2+, Mg2+, Cl−, SO42−, and HCO3−. Due to intense evaporation, salt precipitation can occur in the groundwater; therefore, the possible mineral facies include calcite (CaCO3), dolomite (CaMgCO3), magnesite (MgCO3), gypsum (CaSO4), rock salt (NaCl), and potash (KCl) as shown in Table 2.
3.2. Governing Equations
The basic governing equations for solving the multicomponent motions (mass and energy conservation equations) can be written in the following general form [30]:
The integration is over an arbitrary subdomain of the flow system under study, which is bounded by the closed surface . The quantity in the accumulation term (left-hand side) represents the mass or energy per volume, labels the mass components (e.g., water, air, H2, and solutes), denotes the heat “component,” denotes the mass or heat flux, denotes sinks and sources, and is a normal vector of the surface element that points into .
The simulations were carried out using the reactive transport code TOUGHREACT [23], which introduces reactive geochemistry into the multiphase fluid and heat flow code TOUGH2 V2 [30]. TOUGHREACT is a thermal-physical-chemical code applicable to one-, two-, or three-dimensional geologic systems with physical and chemical heterogeneity. The numerical method for the fluid flow and chemical transport simulation is based on the integral finite difference (IFD) method for space discretization. The system of chemical reaction equations is solved on a grid-block basis by Newton-Raphson iteration. The thermodynamic data used in the simulations were taken from the EQ3/6 database, which was derived using SUPCRT92 [31]. The local equilibrium constants and kinetic rates used in TOUGHREACT are given in Xu et al. [23].
The change in volume of the high concentration solution is related to many factors. Currently, there are three main methods to calculate the density of natural water based on the concentration [32–34]. The relationship between the variations in density and TDS under natural water evaporation is calculated using three methods as shown in Figure 3. There is a relatively significant difference between the results calculated using the third method and the results of the other two methods; in particular, the difference increases gradually when the TDS is greater than 50 g/L. The densities calculated using the first and second methods form essentially linear relationships with the TDS. However, the densities calculated using the third method indicate a nonlinear relationship with TDS. Because the third method calculates the density based on the partial molar volume of the solute and the results calculated using this method are closer to the actual values, they are more suitable for calculating the density of high-concentration brine [34]. In this study, the third method is used to calculate the density of fluids with high TDS concentrations.
3.3. Spatial and Temporal Discretization
We adopt the gradient mesh size method to divide the simulation domain. The mesh size increases gradually from the ground surface to the bottom. The thickness of the smallest cell grid is 0.1 m. The model is divided into 3150 computational units (Figure 4). The time discretization uses the automatic time step method. The minimum time step is 1 s, and the maximum time step is 1000 a. The numerical model will automatically set the time step based on the iterative situation.
3.4. Initial and Boundary Conditions
River infiltration is set as the only groundwater recharge. It is introduced into the model in the form of a specified flow. The infiltration at the piedmont is relatively high, and the infiltration in the lower reach gradually decreases. The spring is calculated using the DELV module, which is similar to the DRAIN module of MODFLOW. The EVAP module is used to simulate groundwater evaporation, and the value of groundwater evaporation is calculated based on the soil water content method [35]. The pressure distribution of the groundwater under stable conditions is taken as the initial pressure condition for the numerical model.
The solute concentrations of well CK7 are taken as the initial chemical conditions for the numerical model (Table 3). The sampling and chemical analysis of the groundwater in well CK7 was conducted in 2013. Its solute concentrations are similar to those of the other wells in the recharge region as well as the river water, even though the ages of the groundwater in these wells vary from thousands of years to 10,000 years. These observations indicate that the chemical components of the groundwater in the recharge region have not changed significantly over geologic time [25]. Therefore, it is reasonable to take the solute concentrations of well CK7 as the initial chemical conditions. Because the groundwater recharge is mainly from river water, the solutes from the river water are taken as the flux boundary conditions. In the discharge region, groundwater discharges by evaporation, but the solutes are left in the aquifer. Therefore, the chemical boundary condition in the discharge region is defined as a no-flux boundary.
4. Results and Discussion
4.1. Validation of the Simulation Model
The simulation results were validated by comparing the simulated groundwater levels with the observed groundwater levels and the simulated solute concentrations with the measured solute concentrations in the groundwater. The observed groundwater levels are from some wells including CK7 and CK8. And the measured solute concentrations are from well CK8. A comparison of the simulated and observed groundwater levels shows that they are similar (Figure 5). The simulated mean residual is 5 m, and the absolute residual is 5.7 m. The regions with relatively large errors are mainly located in the recharge area, where the variation in the groundwater levels is much larger than that in the discharge area. Figure 6 shows a comparison of the simulated TDS with the measured TDS in well CK8. Although the variation in the simulated TDS is smaller, they have the same trends from the shallow to the deep aquifers. The TDS of the groundwater decreases from the surface to a depth of 100 m and then increases. Overall, the simulation model of the groundwater flow and solute transport is reasonable.
4.2. Groundwater Flow System Distribution
Figure 7 shows the simulated distribution of the flow velocity under the present conditions. There are considerable differences in both the horizontal and vertical seepage velocities. In general, the velocity at the top of the basin is greater than that at the bottom of the basin, and the velocity at the basin margin is greater than that in the centre. The maximum difference in the flow velocity can reach several orders of magnitude. The shallow velocity at the margin of the basin is greater than 1 × 10−5 m/s. Groundwater is recharged at the piedmont and is then discharged in the form of springs or surface evaporation at the spill belt. The renewal rate of that groundwater is relatively fast, and it can be defined as the local groundwater flow system shown in Figure 7 (I1, I2). The spring grows into a river, and a portion of the river water infiltrates into the groundwater to form a local flow system (I3). The flow velocities below the local flow system are generally in the range of 1 × 10−6 m/s to 1 × 10−7 m/s, and they are significantly lower than those in the local flow system. Moreover, the seepage depth and seepage distance both increase significantly, although the groundwater still receives recharge from the piedmont. However, this water flows below two local water flow systems (I1, I3) and is discharged via evaporation near the central basin. There are relatively large differences in the flow path from those of the local water flow systems; therefore, it can be defined as the intermediate water current systems (II1, II2). The seepage velocity of the deep groundwater in the central basin is less than 1 × 10−7 m/s, and the flow velocity is extremely slow. The groundwater is supplied from the piedmont and infiltrates to the bottom of the basin. It flows below the intermediate flow system and is discharged via evaporation in the central basin. In comparison to the intermediate flow system, the seepage depth is greater, and the seepage path is longer, so it can be defined as the regional flow systems (III1, III2). In summary, the simulation results indicate that the groundwater in the Qaidam Basin under the present conditions forms a three-level groundwater flow system.
The distribution of these flow systems is controlled by the form of the basin, the sediment permeabilities, and the recharge and discharge conditions. The elevations decrease gradually from the piedmont to the central basin, and the main source for the Qaidam Basin is the river in the proluvial fan. These two factors cause the groundwater to accumulate more potential energy in the piedmont and drive the groundwater towards the lower reaches. The piedmont terrain on the southern side of the basin is higher than that on the northern side, and the recharge is also higher, which causes the development region of the water flow system on the southern side to be larger than that on the northern side. The permeability in the central basin is poor, which causes the groundwater to not be completely discharged from the central basin; instead, it is partially discharged in the form of an overflow spring in the spill belt between the piedmont and the central basin. The river formed by the overflow spring infiltrates again and forms the local water flow system, and it will eventually develop a multilevel flow system.
4.3. Solute Accumulation and Evolution
The simulation results indicate that the groundwater in the eastern region of the Qaidam Basin has been continually undergoing concentration and enrichment since 70 ka BP (Figure 8). However, with climate change, the distribution of the groundwater flow system and the solute distribution may change.
The climate conditions at 70–45 ka BP were similar to the current conditions. Due to the influence of evaporation, the hydrochemical components were first concentrated in the local flow system at the spill belt. This mainly occurred because of the large water circulation in the spill belt, where the water supply was adequate and the evaporation was high; therefore, the concentration first occurred in the spill belt, which resulted in salty and brackish water. However, the hydrochemical components concentrated in the spill belt did not remain in the spill belt. They were transported with the water towards the central basin and were eventually enriched in the central basin. Later, at 30 ka BP, the groundwater components continued to be enriched. By 45 ka BP, most of the groundwater in the central basin had evolved to saline water to form a distribution with salty water in the upper part of the overflow zone and freshwater in the lower part.
The time period of 45–25 ka BP was a minor interglacial period; it was the relatively warm time of the last glacial period. More rainfall fell than at present, and the potential evaporation was relatively low [26]. The reduction in evaporation resulted in greater groundwater discharge in the form of springs, and the local flow system was enhanced. The discharge centre of groundwater shifted towards Nuomuhong, and the enrichment centre of the water-soluble components also moved. The relative changes in the circulation of the different flow systems also caused the distributions of the shallow freshwater and deep salty water to expand; namely, some desalination of the shallow water occurred.
The time period of 25–15 ka BP was the Pleniglacial during the last glacial period, and the rainfall was only 30% of the present level. The abrupt reduction in recharge caused the enhancement of the regional water flow system, and the water and salt enrichment centres returned to the central basin. The variation in the density difference between the shallow water and deep water generated free convection of the groundwater in the central basin, and the groundwater with relatively low concentrations moved upward via buoyancy, whereas the groundwater with relatively high concentrations moved downward. However, the free convection caused by this density difference was only temporary over geologic time, and the forced convection driven by gravity rapidly became a dominant factor.
The period from 15 ka BP to the present includes the late glacial period and postglacial period. The temperature gradually increased again, and the recharge and potential evaporation gradually changed to their present levels. The distributions of brackish and freshwaters were similar to those in the initial stage of the simulation. The distribution of saline water reached its maximum and formed a wedge of freshwater below the spill belt.
Currently, the groundwater in the central basin is mostly saline water with a concentration of 10–20 g/L. If there are no significant changes in the climate conditions in the future, the soluble components of the groundwater will continue to be concentrated and enriched. The calculated results of mineral saturation indicate that under the current conditions, except for dolomite, most minerals have not reached the saturated state. Currently, the groundwater is in the carbonate precipitation stage, and it is still far from precipitating gypsum and halite.
Due to the enrichment of salt in the central basin, the area of high TDS grew over time (Figure 9). Initially, the area of increased dramatically, which was caused by groundwater evaporation. The salt was transported from the surface of the central basin towards the bottom of the basin and from the centre of the basin towards the margin. After 65 ka BP, when the groundwater with high TDS occupied the largest area of the central basin, the area of continued to increase but at a relatively slow rate. This occurred because the salt transport from the central basin (discharge region) towards the margin of the basin (recharge region) was dominated by hydraulic dispersion; however, when the salt was transported close to the recharge region, hydraulic dispersion was prevented by the groundwater flow from the recharge region. The areas of higher TDS, such as and , increased much more slowly than the area of . Over time, larger areas became occupied by groundwater with higher TDS. Currently, the TDS in 81.7% of the area is higher than 3 g/L, and 55.0% of the area is higher than 10 g/L.
4.4. Variation in Chemical Composition
There are significant variations in the hydrochemical components and composition from the recharge region to the discharge region (Figure 10). The variations in the concentrations of different components were not large before 45 km and after 70 km (between the alluvial fan and the spill belt), which is consistent with the initial supply concentration. The main cation is Na+, and the main anion is Cl−; this is mainly because the chemical action in this region is dominated by leaching and mixing. The sediment has experienced sufficient filtration over geologic time; therefore, the hydrochemical components will no longer significantly increase. The components in the groundwater mainly come from surface water recharge. Between 45 km and 70 km (in the central basin), except for the HCO3− concentration, which decreases, the concentrations of all the other components increase significantly. However, there are still differences in the rates of increase of the different components. The concentrations of Cl− and Na+ increase the most, and their maximum concentrations are 84 and 79 times higher than their initial concentrations, respectively. The concentrations of Ca2+ and Mg2+ increase moderately, and their maximum concentrations are 43 and 51 times higher than their initial concentrations, respectively. The evaporation and concentration mainly occur in the central basin. HCO3− reacts with Ca2+ and Mg2+ to generate dolomite, which causes a decrease in the HCO3− concentration. The formation of dolomite will generate H+, which increases the concentration of H+ in the groundwater and eventually causes the pH of groundwater in the central basin to decrease from 7.7 to 7.2. In general, evaporation enhances the concentrations of most components, but there are no significant changes in the proportions of the different components. This occurs mainly because the concentration has not resulted in a period of significant diagenesis.
5. Conclusions
Using the reactive transport simulation technique, a vertical 2D numerical model was developed to simulate groundwater flow and solute transport since the last glacial period. The simulation was performed on a typical cross section across the Qaidam Basin, which is located in the northeastern Tibetan Plateau.
The Qaidam Basin generally receives recharge from the piedmont regions and is drained by evaporation. Due to the basin morphology, sediment lithology, and recharge and discharge conditions, a three-level nested groundwater flow system is developed. There are significant differences in the flow ranges and velocities of the different groundwater flow systems. In the local water flow system, the seepage path and seepage depth are relatively short, and the flow velocity is relatively fast. In the regional water flow system, the seepage path is relatively long and the seepage is relatively deep (it can extend to the bottom of the basin). The flow velocity in the regional is relatively slow, which can be many orders of magnitude slower than that in the local water flow system.
Since the last glacial period, the groundwater in the eastern part of the Qaidam Basin has experienced solute concentration and enrichment. In general, the solute is recharged with the groundwater in the piedmont regions and moves towards the central part of the basin. Due to intense evaporation and concentration, the solutes are concentrated and enriched in the central basin. During different geological periods, climate change affected the distribution of the groundwater flow system and the solute concentrations. During warm climate periods, the local flow system was enhanced, and the drainage centre of the groundwater and the solute concentration centre both shifted towards the side with the greater recharge. During periods of drought, the drainage centre and solute concentration centre shifted towards the central part of the basin. Currently, the groundwater is under conditions favourable for the precipitation of carbonate but is still far from the precipitation of gypsum or halite.
Many parameters are involved in the coupled model of water flow and solute transport. In particular, for a large basin such as the Qaidam Basin, the hydrogeological parameters of the deep aquifer are not easy to obtain, and the hydrodynamic dispersion coefficient under large-scale conditions is also relatively difficult to determine. Therefore, further analysis of the numerical model parameters for the water flow and salt migration is necessary. In addition, the study did not take into account the influence of human activities on the groundwater flow and salt migration, which also need to be explored further.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This work was jointly supported by the Natural Science Foundation of China (no. 41702282), the Fundamental Research Funds for Central Public Research Institutes (no. YYWF201626), and a China Geological Survey Project (no. DD20160238).