ABSTRACT
The Tibetan Plateau (TP) plays a vital role in supplying water to approximately two billion individuals residing downstream, as it serves as the source of numerous significant rivers across Asia. Nevertheless, the significance of groundwater in sustaining the plateau rivers remains uncertain, hence jeopardizing the long-term water resource sustainability of Asia. Here, we assess the groundwater contribution to river flow in the middle reaches of the Yarlungzangbo River (YR), a typical TP location, and we pinpoint the vital role that plateau groundwater plays in sustaining the water balance of the river. The findings indicate that the YR mainstream has a runoff excess of approximately 36–42 × 108 m3/a, which accounts for 22.8–26.5% of the mainstream flow and is strongly correlated with groundwater discharge. The research findings also indicate that the movement of groundwater on the plateau is primarily influenced by the variation in land elevation and replenishment at higher elevations. The contribution of groundwater is inversely proportional to the topographic steepness and precipitation, and its recharge mechanism is associated with the large-scale active fault zones in the vicinity of the YR.
HIGHLIGHTS
Groundwater supports the base flow of the Yarlung Zangbo River (YR), providing 22.8–26.5% annually to the middle reaches of the Yarlungzangbo mainstream.
The groundwater flow in the plateau is driven and sustained by the topographic gradient and recharge at high elevations.
A series of tensile fault structures in the YR basin are potential channels for groundwater circulation and discharge.
INTRODUCTION
The Tibetan Plateau (TP) serves as the origin of numerous major rivers across Asia and plays a crucial role in providing freshwater to the extensive low-lying regions of China, Bhutan, Nepal, and India (Bookhagen 2012). The elevated terrain of the TP functions as a crucial source of freshwater for a significant portion of the global population. However, climate change is threatening this massive reservoir (Huss & Hock 2018). Due to an increase in rainfall and a decrease in frozen precipitation, the water levels of various lakes, rivers, and underground water sources on the TP and in other mountainous regions have increased significantly since the late 1990s (Langman et al. 2022; Yao et al. 2022; Lorenzi et al. 2024). The groundwater reserves grew at a rate of 5.0 × 109 m3/a from 2003 to 2009, a growth rate that is similar to that of lake water storage (Jiao et al. 2015). Hence, it is of utmost importance to understand the mechanisms through which groundwater originates and its influence on rivers across the TP.
The groundwater system has a crucial function in sustaining river flow during the dry season and storing precipitation throughout the monsoon season (from late May to September), thereby mitigating the peak river flow (Cuthbert et al. 2019). The proportion of groundwater in alpine rivers surpasses that of snow and glacial meltwater (Bookhagen 2012). For instance, groundwater plays a crucial role in sustaining the flow of rivers downstream in the Pamir Mountains, particularly during the winter season when groundwater serves as the primary water source. Groundwater accounts for approximately 40% of the total annual river flow (Pohl et al. 2015). However, the water budget of TP rivers is usually studied in terms of surface water composition, and the consensus is that precipitation and the thawing of ice and snow are the primary contributors to maintaining water equilibrium. The influence of groundwater on rivers is seldom considered, greatly limiting the precision of water resource evaluations within the basin (Li et al. 2020). While the contribution of groundwater discharge to rivers in terms of available water resources is acknowledged, there remains limited understanding regarding the extent to which rivers on the TP rely on groundwater. Recent research has investigated the significance of groundwater systems in mountainous regions on the TP, specifically focusing on faulted areas. Additionally, deep circulating groundwater originating from thermal springs can also contribute to river flow upon emerging at the surface (Andermann et al. 2012). This particular subsurface water movement is typically propelled and sustained by variations in land elevation and replenishment in elevated areas through the infiltration of precipitation and melted snow (Ge et al. 2008; Schmidt et al. 2020). However, there is still a lack of comprehensive knowledge regarding the mechanisms through which groundwater flows into rivers. The quantification of the recharge discharge process and the contribution of groundwater to river discharge remain insufficiently determined. These essential yet absent details are crucial for predicting the future availability of water and alterations in the ecohydrological system of the TP.
The Yarlung Zangbo River (YR) is a representative example of a major river on the TP, and it flows along the collisional boundary between the Indian Plate and the Eurasian Plate (Shen et al. 2022). The most extensive transverse rift zone in the basin is located in the middle reaches of the YR. Significant tectonic activity has the potential to modify the flow pattern of groundwater controlled by topography, leading to an increase in the amount of groundwater contributing to runoff. Consequently, this alteration may impact the composition of rivers. However, previous relevant studies focused mainly on the surface water cycle level, and the abnormal variation in river isotopic compositions was thought to be caused by the spatial zoning of precipitation-forming water vapor and that this was only marginally influenced by the considerable contribution of groundwater to runoff (Yao et al. 2022). Alternatively, other studies have focused on the vertical circulation level of groundwater and have found that groundwater aquifers are recharged by atmospheric precipitation infiltration and discharged by evaporation, concluding that the discharge of shallow groundwater causes changes in the isotopic composition of rivers (Ge et al. 2008). The YR mainstream exhibits an atypical composition during winter despite the dissimilarity in precipitation and water vapor, which can be attributed to more than 90% of the basin's rainfall occurring during the summer season. In alpine basins, the primary process responsible for recharging and discharging groundwater is the vertical and horizontal circulation of water within a specific structure. Therefore, the patterns of groundwater flow under tectonic constraints, as well as its impact on rivers, remain poorly constrained. Our previous research revealed anomalous physical and chemical characteristics in the isotopes, hydrochemistry, and water temperature of the middle reaches of the YR. It was concluded that high-elevation surface water enters the groundwater system through fault zones and is discharged into the YR valley. However, to date, the contribution of groundwater to the runoff of the YR remains undocumented, and the dynamic mechanisms governing groundwater recharge are not yet fully understood.
Therefore, we aimed to examine the hydrological mechanisms of plateau groundwater by focusing on the middle reaches of the YR as a case study. The systematic use of physical and chemical factors, Rn isotopes, and numerical simulation methods reveals the coevolution mechanism of geomorphic structure and groundwater recharge–runoff–discharge processes in alpine valleys, with the goals of (i) identifying the groundwater flow patterns and spatial distribution characteristics of typical alpine basins and (ii) reevaluating the groundwater discharge and its contribution to the runoff of the middle YR, (iii) uncovering the mechanism behind groundwater circulation in typical alpine basins. The originality of this paper lies in the coupling of isotopic and numerical modeling methodologies in order to quantify the crucial role that plateau groundwater plays in sustaining river base flow.
GEOGRAPHICAL AND GEOLOGICAL SETTING
(a) Geographic location of the TP and YR basin. (b) The YR basin's terrain is depicted with a color scheme that corresponds to elevation using a 1 km digital elevation model. The labels indicate the positions of the sample sites employed in this research. (c) Estimated distribution of various categories of permafrost and seasonally frozen ground in the study area.
(a) Geographic location of the TP and YR basin. (b) The YR basin's terrain is depicted with a color scheme that corresponds to elevation using a 1 km digital elevation model. The labels indicate the positions of the sample sites employed in this research. (c) Estimated distribution of various categories of permafrost and seasonally frozen ground in the study area.
According to the formation lithology and aquifer media, the groundwater in the study area is primarily quaternary loose rocks pore water and bedrock fissure water (Liu 2014). The quaternary loose rocks pore water is primarily distributed in the valley and fore-mountain zone in the intermountain region, with the lithology consisting primarily of gravel and gravel layers with loose structures. The depth of the groundwater level is generally <30 m, the aquifer thickness ranges from several meters to hundreds of meters, the hydraulic conductivity is 10−6–10−3 m/s, the specific storage is (4.9–10) × 10−4 m−1, and the specific yield is 0.025–0.35. The bedrock fissure water is primarily found in the bedrock mountain ranges south and north of YR. The aquifer is made up of granite, biotite monzonitic granite, and diorite granite, which are medium acid and ultrabasic intrusion rocks that range from the Yanshanian to the alpine stage. The rock is relatively broken, the hydraulic conductivity is 10−8–10−11 m/s, the specific storage is <3.3 × 10−6 m−1, and the specific yield is 0.001. The rock mass topography is steep, with a thin weathering crust. Furthermore, ground vegetation is sparse, topographic relief is extensive, precipitation is excreted smoothly, and residence duration is brief. Therefore, groundwater recharging takes a short period.
MATERIALS AND METHODS
Samples and analyses
In May 2018, field observations and sampling activities were carried out. The electrical conductivity (± 1.0%), pH (± 0.1), and temperature (± 0.1°C) were measured in the field at the time of sampling using calibrated hand-held field meters (HACH: HQ40D). 222Rn was analyzed in the field using a DURRIDGE USA (RAD 7) radon detector. The reported data were averaged from five measurements taken within approximately 1 hour. The detection limit was 0.4000 cpm/pCi/L, and the analytical precision was within ±5% (1σ). The sampling points are shown in Figure 1(b).
Model data
The 90 m Shuttle Radar Topography Mission digital elevation model (DEM) was employed to delineate the terrain features and derive the catchment boundaries for the YR basin. The study area shares the same north‒south boundary as the YR basin, while its western boundary is defined by a river, and its eastern boundary is demarcated by another subbasin (Yangcun Basin). The high topography between these two subbasins serves as a natural hydrological divide (Figure S1). The monthly precipitation and temperature records for Lazi, Rikaze, and Nugesha from 1982 to 2022 were acquired from the China Meteorological Information Center (http://data.cma.cn/). The monthly gridded precipitation product of the TP (0.1° × 0.1°, http://data.tpdc.ac.cn) served as a means of validating the interpolated precipitation of the YR Basin, as the meteorological stations are all located below 4,500 m above sea level (Figure S2). The monthly runoff data from the Rikaze and Nugesha stations for the period of 2013–2018 were acquired from the Hydrology and Water Resources Survey Bureau of Rikaze. Daily groundwater level data for 12 boreholes from 2018 to 2022 were obtained from the China Institute of Geological Environment Monitoring. Data on the distribution of glaciers and permafrost in China were obtained from the second dataset of China's Glacier Inventory.
Groundwater model
The numerical groundwater model is constructed based on the Finite Element subsurface FLOW system (FeFlow) software developed by WASY, a German company. The finite element method is employed to discretize the structure, which transforms a complex object into a finite number of easily analyzable units. Additionally, the unknown function within the study area is approximated by employing fragment functions on each unit. The unit discretization process was applied to the model, resulting in the division of the study area into 2,510,236 triangular units. These units were further characterized by a total of 1,272,392 nodes (Figure S3). The vertical dispersion is five layers with a total thickness of approximately 3.5 km. The initial model was developed in a state of equilibrium to establish the fundamental pattern of subsurface water movement within the research region. Recurring models were subsequently created using daily time intervals, building upon the established steady-state model to evaluate the yearly and seasonal proportions of groundwater discharge concerning overall river currents (detailed information on model development is provided in the Supplementary Material).
The model enables the assignment of time series parameters to accurately simulate the changeable climate and hydrological conditions in reality. Therefore, time series parameters are used when parameter assignments are made. In winter (November–April), the seasonally frozen aquifer can be regarded as a weakly permeable layer (Figure 1(c)), with an average hydraulic conductivity (AHC) of 8.2 × 10−10 m/s. The main component of the permafrost unit was envisioned as a consistently weakly permeable layer with an AHC of 5.2 × 10−11 m/s (Yao et al. 2021). The alpine permafrost layer was envisioned as a somewhat less cohesive, weakly permeable layer than the northern uninterrupted permafrost layer, possessing an AHC of 6.5 × 10−11 m/s (Figure S4). In summer (May–October), the seasonally frozen layer melts, forming an aquifer with an AHC of 3.0 × 10−4 m/s, while the thawed areas outside the river valleys exhibit a lower AHC of 7.7 × 10−6 m/s. The ‘Discrete Feature’ module was used to generalize the fracture zone with an AHC of 10−3 m/s. Darcy's law was applied to define the fluid movement within the fault. Groundwater recharge was determined as the product of precipitation and the rainfall infiltration coefficient (α) for different months. The precipitation infiltration coefficient is based on the empirical values of different land distribution types. The α values are 0.4–0.8 in the river valley area, 0.1–0.2 in most areas of the mountain area, and 0.2–0.4 in the alpine source area, which is characterized by extensive wind compounds and moraine-covered regions with significant permeability. The initial values of recharge and hydraulic conductivity in the study area were obtained from Ge et al. (2008), Chen et al. (2021), and Yao et al. (2021) and adjusted by manual calibration and the PEST inversion algorithm (Doherty & Hunt 2010). Relevant data from the study area, including winter base flow records (from November to April) and borehole groundwater level monitoring data, were utilized as calibration targets for the transient model.
RESULTS AND DISCUSSION
Model calibration
River flow primarily originates from precipitation, glacial ice melting, and the discharge of groundwater (Andermann et al. 2012). However, the majority of river flow is attributed to glacial meltwater during the summer (May–October), while its contribution decreases significantly in the dry winter period (November–April). As glaciers are situated in both the southern and northern parts of the basin and are located in the upper source region, the significance of glacial meltwater is limited to this specific area. Furthermore, the amount of winter precipitation primarily consists of snowfall and is less than 10 mm per month. Therefore, the recharge of river flow by surface runoff derived from precipitation and meltwater is assumed to be negligible in the dry winter months. During the winter season, the main contributor to river flow in the middle reaches of the YR basin is groundwater (Hayashi 2020). Stream baseflow, derived from groundwater discharge, can be most accurately estimated by analyzing the observed streamflow data during this period (Walvoord et al. 2012). Therefore, the model parameters are adjusted and calibrated based on the base flow in the Nugesha and Rikaze subbasins in winter and the groundwater level monitoring data at each borehole point.
Simulated and observed (during 2018–2022) borehole water levels for the steady-state model (the location of drilling points is shown in Figure S5).
Simulated and observed (during 2018–2022) borehole water levels for the steady-state model (the location of drilling points is shown in Figure S5).
Contrast between simulated and observed data. (a) Water levels in pumping wells within the study area were compared between computed and observed values under steady-state conditions (the position of the observation point is shown in Figure S5). (b) Baseflow values were determined and compared between the two hydrological stations in relation to the steady-state model. The blue triangles represent the baseflow calculated by the dynamic model, while the boxplot displays the observed streamflow during periods of low winter precipitation.
Contrast between simulated and observed data. (a) Water levels in pumping wells within the study area were compared between computed and observed values under steady-state conditions (the position of the observation point is shown in Figure S5). (b) Baseflow values were determined and compared between the two hydrological stations in relation to the steady-state model. The blue triangles represent the baseflow calculated by the dynamic model, while the boxplot displays the observed streamflow during periods of low winter precipitation.
Groundwater flow patterns
(a) Geographical distribution of groundwater flow. (b) The simulated average annual baseflow discharge (Qg) per kilometer along streams and groundwater recharge (R) per unit basin area.
(a) Geographical distribution of groundwater flow. (b) The simulated average annual baseflow discharge (Qg) per kilometer along streams and groundwater recharge (R) per unit basin area.
Groundwater recharge is the process by which an underground water-bearing system receives water from the surrounding environment. This term refers to the recharging of underground water via precipitation, meltwater, irrigation, underground runoff, or river seepage. Groundwater discharge refers to the release of water from subterranean aquifers into the surrounding ecosystem, and springs serve as natural locations where this discharge occurs. Hot springs are widespread in the middle reaches of the YR (Figure S1), and groundwater discharge varies greatly from less than 0.1 to more than 30 × 106 m3/km/a (expressed in unit river length). The maximum discharge is mainly concentrated in the river valley where the north–south-trending fault zones meet the YR, reflecting the water-transport effect of the north–south-trending tensile fault zones. The study area experienced an increase in replenishment, with the highest rates observed in the central YR valley and lower rates observed toward the south and north, from 0 to 450 mm/a (expressed in unit basin area). In the middle of the YR, in 82.3% of the area, the groundwater discharge is less than 10 × 106 m3/km/a, and in 78.2% of the area, the groundwater system recharge is less than 200 mm/a (Figure 6((b)). Areas with high groundwater discharge to rivers (>10 × 106 m3/km/a) are concentrated in regions characterized by elevated recharge rates (300–500 mm/a). These areas are primarily situated within the alluvial sector at the confluence of the Rikaze and Nugesha subbasins. The significance of this location stems from its role as a central hub for population and agriculture in the middle reaches of the YR. Consequently, irrigation practices and reservoir management contribute significantly to augmenting both recharge and discharge rates. Furthermore, based on the simulation results, the YR mainstream area from the Lazi to Nugesha stations has a groundwater discharge of 15–20 × 106 m3/km/a. The river reach is approximately 200 km long, meaning that approximately 30–40 × 108 m3 of groundwater is discharged annually to sustain the base flow of the YR. This result confirms the significant contribution of groundwater to the flow of the middle YR, as it is similar to the groundwater discharge calculated in our earlier study (36 × 108 m3). Further discussion regarding the water source and the circulation mechanism of such a large discharge is still required.
Implication of physicochemical parameters and FBM and 222Rn MBM
(a) Spatial variation in water temperature in the YR middle reaches; (b) spatial variation in pH in the middle reaches of the YR.
(a) Spatial variation in water temperature in the YR middle reaches; (b) spatial variation in pH in the middle reaches of the YR.
The calculation parameters and results of the FBM and 222Rn MBM
FBM . | Longitude/(°E) . | Latitude/(°N) . | Runoff/Q (m3/s) . | Precipitation/P0 (mm) . | catchment area/S (×104 km2) . | Groundwater discharge/Vg (×108 m3) . |
---|---|---|---|---|---|---|
LZ | 87.68 | 29.16 | 151.46 | 370.7 | 6.720 | 35.982 |
NGS | 89.78 | 29.38 | 500.27 | 435.6 | 13.52 | |
MBM . | 222Rn concentration/C (Bq/m3) . | Flow velocity/v (m/s) . | Depth/h (m) . | Length of reach/l (m) . | groundwater discharge per unit length/q (×10−4 m3/s/m) . | Groundwater discharge/Vg (×108 m3) . |
LZ | 485.07 | 4.5 | 4 | 340.19 | 3.902 | 41.865 |
NGS | 708.83 |
FBM . | Longitude/(°E) . | Latitude/(°N) . | Runoff/Q (m3/s) . | Precipitation/P0 (mm) . | catchment area/S (×104 km2) . | Groundwater discharge/Vg (×108 m3) . |
---|---|---|---|---|---|---|
LZ | 87.68 | 29.16 | 151.46 | 370.7 | 6.720 | 35.982 |
NGS | 89.78 | 29.38 | 500.27 | 435.6 | 13.52 | |
MBM . | 222Rn concentration/C (Bq/m3) . | Flow velocity/v (m/s) . | Depth/h (m) . | Length of reach/l (m) . | groundwater discharge per unit length/q (×10−4 m3/s/m) . | Groundwater discharge/Vg (×108 m3) . |
LZ | 485.07 | 4.5 | 4 | 340.19 | 3.902 | 41.865 |
NGS | 708.83 |
Note. LN refers to Lazi station; NGS refers to Nugesha station.
(a) Spatial variation in runoff in the YR; (b) spatial variation in Rn isotopes in the middle reaches of the YR.
(a) Spatial variation in runoff in the YR; (b) spatial variation in Rn isotopes in the middle reaches of the YR.
Control of groundwater flow pattern and discharge
Relationships between the maximum variations in elevation along flow paths and their corresponding distances.
Relationships between the maximum variations in elevation along flow paths and their corresponding distances.
A 3D model of the YR where the color palette represents elevation. (a) Groundwater level in the middle reaches of the research site; (b) conceptual model of groundwater circulation in the study area.
A 3D model of the YR where the color palette represents elevation. (a) Groundwater level in the middle reaches of the research site; (b) conceptual model of groundwater circulation in the study area.
In addition to topographic control, basin fault structure has a considerable impact on groundwater flow. The groundwater discharge in the middle reaches of the YR, particularly in the valley area, is greater than that of large rivers on plains. Over a river length of less than 200 km, approximately 36–42 × 108 m3 of groundwater is discharged annually to maintain the base flow of the YR. If such a large flow is realized only by the slow infiltration of precipitation along the unsaturated zone in traditional hydrology, the whole cycle would take a long time, which contradicts the groundwater age and the weak water‒rock interaction of the groundwater observed in a previous study (Tan et al. 2014). Therefore, the groundwater replenishment of the YR mainstream may include not only water moving slowly through the unsaturated zone or surface fractures but also concentrated groundwater drainage with rapid circulation of alpine meltwater and precipitation along deep and large fractures in the middle reaches of the YR. Recent research suggests that deep groundwater may represent the world's largest ‘reservoir’. Because the majority of this groundwater is deep and is stored in rocks with low permeability, it has difficulty circulating or flowing to the surface. If the high-elevation recharge areas are proximal to the low-elevation discharge areas and a significant fault is present, substantial pressure differentials will occur. Consequently, this pressure disparity will facilitate the movement of shallow water toward deeper regions, where it can actively participate in deep circulation before eventually resurfacing through water conduction within the fault zone (McIntosh & Ferguson 2021; Shi et al. 2021). Extensive deep faults with nearly north‒south orientations were observed in the study region, particularly in the southern section of the basin (Figures 2 and S1). In terms of mechanical characteristics, these fault structures are mostly defined by tension or tensional shear, which considerably aids the circulation of precipitation and meltwater from the southern Himalayas to the YR valley via fault zones (Figure 10(b)). The distance of lateral migration of high-elevation recharge water along a fault may exceed 150 km (Su et al. 2020), and the water flow in each extensional fault will be driven by the large difference in the water head generated by the topography and by the density difference between hot and cold water. The recharge water will be discharged as spring water at the fault's lowest point where the fault meets the YR.
CONCLUSIONS
An FeFlow-based groundwater flow model was developed to investigate the dynamics of groundwater flow in the Lazi–Nugesha reaches of the YR basin. The model considers a variety of elements, including precipitation, fault structure, and permafrost, and combines river water physicochemical properties and Rn isotopes to identify groundwater discharge mechanisms. The movement of groundwater in the examined region is related to variations in topography, in which a greater disparity in elevation leads to a reduced distance traveled by the flow. In the middle reaches of the YR, short (<10 km) flow paths account for more than 70% of the spatial flow pattern of the groundwater system, indicating rapid groundwater renewal. The recharge of groundwater occurs mainly through precipitation or meltwater at high elevations, followed by discharge in valleys and faults. Furthermore, data on water temperature, pH, and Rn isotopes indicate the driving effect of water conduction in fault zones on groundwater recharge and discharge. The annual discharge of the groundwater that reaches the mainstream of the YR along fractures is approximately 36–42 × 108 m3, accounting for 22.8–26.5% of the discharge at the Nugesha station. The significance of groundwater in sustaining the continuous flow of alpine rivers is underscored, facilitating investigations into regional hydrological phenomena. However, due to the complexity of actual samples and the impact of simulation parameters, the estimation of groundwater's contribution to YR flow may be biased. Therefore, additional research into groundwater migration patterns is required in the future to provide more accurate groundwater resource assessments and planning.
AUTHOR CONTRIBUTIONS
D. S. conceptualized the process, designed the work, contributed in the data collection, and wrote the original draft, developed the methodology, and provided the software; S. H. wrote the review and edited the article. All authors have read and agreed to the published version of the manuscript.
FUNDING
This study was financially supported by the School Scientific Research Development Fund project (2022LFR091).
ACKNOWLEDGEMENTS
We appreciate Yu Zhang, and Peixin Cong for their assistance with the field work. We also thank the editors for their valuable comments.
DATA AVAILABILITY STATEMENT
All relevant data are included in the paper or its Supplementary Information.
CONFLICT OF INTEREST
The authors declare there is no conflict.