With rapid urbanisation, a karst water recharge area of the Jinan spring catchment was damaged. Thus, managed aquifer recharge projects were built in the western Jinan spring catchment to protect the water supply of the spring. Yufuhe River was selected as the study area to compute the effective recharge rate into karst aquifers. This strong seepage zone has a large gradient and undergoes a specific hydrogeological condition in which two strata of a gravel layer and limestone change to three strata of gravel, impermeable clay shale and limestone at the open window of the karst aquifers. A hydraulic model called HEC-RAS was applied to simulate the river stage, and a numerical groundwater model called HYDRUS-3D was adopted to simulate the groundwater mound dynamics and estimate river flow seepage into the aquifers. The effective recharge rates are 64.9%, 65.2% and 68.1% when the buried depths of groundwater are 40, 30 and 25 m. An analysis of the electric conductivity, water table, temperature and water volume data found an effective recharge rate of 68.3%. Results of field monitoring confirmed the accuracy of the numerical simulation and showed that most of the recharged water in the study reach can be effectively recharged into the karst aquifers.
Managed aquifer recharge (MAR) refers to the purposeful recharging of water in aquifers to accelerate their recovery or to benefit the environment (Dillon et al. 2009). Wang et al. (2015) treated rainwater on the basis of the groundwater quality standard and recharged it into fracture karst aquifers. The recharged water was used for drinking and protecting the groundwater environment. Water from the Yellow River was recharged into the ground through an underground-channel pipe-shaft system to increase agricultural irrigation water and restore a groundwater over-exploitation funnel (Rong et al. 2016). However, designing and testing an MAR project is time-consuming and costly. Numerical models are suitable for evaluating MAR projects with the effects of the complexity of hydrological systems, such as temporally varied boundary conditions (e.g. evapotranspiration from a channel bank or unsaturated flow beneath a stream) or complex hydrogeology (e.g. heterogeneous or anisotropic aquifer) (Huang et al. 2015). Sallwey et al. (2018) used HYDRUS-2D/3D in designing a small-scale seepage pool and discussed the advantages and disadvantages of modelling in MAR evaluation. Contact area and discharge time of water and river channels considerably impacted the amount of recharged groundwater during MAR projects. Besides, the exchange in surface water and groundwater is influenced by many factors, such as hydrogeological conditions, topography and climate of a recharge site; therefore, the interaction between surface water and groundwater is a complex process (Gao et al. 2018). A surface water and groundwater coupled model through a numerical simulation is an important method for investigating the water seepage process of river surface water. Sophocleous & Perkins (2000) simulated the severe depletion of groundwater in Kansas by establishing a SWAT and MODFLOW coupled model and proposed three practical management applications. In the present study, a strong seepage area of the Yufuhe River between Dongkema and Cuima villages was selected under specific hydrogeological conditions as the study area to describe the details.
MAR is an effective measure for spring restoration. In the upper reaches of the Yufuhe River, two multi-source water supply projects have been constructed. These projects play an important role in improving surface water quality, increasing domestic use of karst water supply in Jinan, maintaining a high water table of karst water in West Jinan and protecting the spring. A high-efficiency and low-cost recharge scheme for the Yufuhe River channel under specific landform and hydrogeological conditions is addressed. In the numerical simulation method, a surface water and groundwater interaction is realised using a semi-coupled HEC-RAS and HYDRUS-3D model. The discharge time and effective recharge of karst water have been calculated under various water depths when the karst water table reaches the bottom of the Quaternary. In monitoring data analysis, the amount of water that effectively replenishes karst water, surface water that enters the downstream river channel and pore water loss in the river gravel layer are obtained on the basis of surface flow measurement and shallow and deep groundwater-level monitoring analysis and calculation.
The study area illustrated in Figure 1 is located in the Yufuhe River Basin in West Jinan, China. It is a semi-humid region with mean annual precipitation of 670 mm, mainly concentrated during the summer season of June–September. The length of the stream from the Dongkema Village to Cuima Village with a strong seepage is of 5,560 m, and the width is approximately 50 m. The gradient of the streambed is 1/500, with significant characteristics of karst landform.
Jinan is located at the north wing of Mount Tai and is generally a monoclinic structure inclined to the north with Palaeozoic strata as the main body distributed from old to new with a dip of 10°. In Figure 2(a), the stream water does not recharge the karst aquifer until it passes through the study area with strong infiltration. The vertical stratifications of the strata in the study reach are the Quaternary sand gravel layer at 20–25 m, the Zhangxia Formation limestone at 140–180 m and the aquitard Xuzhuang Formation. The aquifer has a minimum thickness on the upstream and increases whilst the downstream flows further away because of water erosion. The open window of the Zhangxia karst aquifer disappears beyond the strong infiltration reach. The invalid recharge to Quaternary pores or water flowing through the surface is caused by the large velocity of water in the streambed which is located in the sloping plane of a piedmont. Therefore, the effective recharge rate of karst water is reduced. The upstream of Yufuhe River is Wohushan Reservoir, which is formed by the convergence of three branches and when water needs be released during a dry period depends on the state of Jinan spring flowing. The area of the control watershed is 557 km2, accounting for 67% of the total area of the watershed. The total design capacity of the reservoir is 118.5 × 106 m3.
Three monitoring sections are set up to monitor the seepage of groundwater in the Yufuhe River discharge project depicted in Figure 1. The first section is located near Dongkema Village, and only one Quaternary pore water monitoring well (Q1) is established due to the Quaternary water overlying shale of a confining bed. The second section is located near Cuima Bridge, including a Quaternary pore water monitoring well (Q2) and a Zhangxia Formation karst water monitoring well (S1). The third section is located near the upper right bank of Wenshan Pumping Station, including a Quaternary pore (Q3) and an Ordovician karst water monitoring well (S2). The groundwater stage, temperature and conductivity can be automatically monitored by remote monitoring devices. Moreover, the rates of released water for each event are recorded by a pump station or reservoir management office.
Numerical simulation method and results
The HEC-RAS model was used to compute the hydrograph, and the river stage was used as a transfer variable to the streambed of the HYDRUS-3D model under a constant water head boundary condition. Therefore, a semi-coupled model of the river–groundwater system is established.
The water balance of the river–groundwater system shown in Figure 2(c) involves three parts, namely, the water balance of the river subsystem, the water balance of the groundwater subsystem and the overall water balance of the subsystem. Given the specific hydrogeological conditions at an open window, the water balance process in the study area is complicated. The water outlet is near the upper reaches of the study area. When the surface water moves towards the study area, the river water begins to infiltrate heavily. Part of the water simultaneously flows downstream along the channel given the influence of the steep riverbed. The remaining amount of released water that penetrated below the ground surface is mostly transferred to Quaternary and karst waters. In the latter, part of the Quaternary water further discharges out of the open window whilst the loss element is influenced by the downstream impervious layer and riverbed slope.
The study area is 4.45 km2 with 800 m width and 5,560 m length. The groundwater domain is conceptualised as a heterogeneous and isotropic geological medium. A 3D groundwater flow model, namely, HYDRUS-3D (Šimůnek et al. 2008), was used to simulate groundwater flow and evaluate the effective recharge amount into the aquifer.
In Figure 2(b), BC represents the constant head boundary, and the head is equal to the river depth of the section calculated using HEC-RAS. GH, GF and EF are generalised to the impermeable boundaries; AH and DE are the free drainage boundaries, and AB and CD are the atmospheric boundaries.
In the modelling, the study area was discretised into two layers in the vertical direction. The top layer consists of sand gravel and loam with 20–25 m thickness, thus underlying the Zhangxia Formation with 140–180 m thickness. The mesh is divided into 3,300 non-uniform triangles with an interval of 50 m and no stretching. The grid is encrypted at the bottom of the riverbed. In accordance with the position of the monitoring well S1, an observation point is set at the same location. A spatial variation in terms of the interpolation of observed groundwater tables was used as an initial condition for groundwater modelling.
The model parameters were initially determined in accordance with pumping tests and were further calibrated using observed groundwater tables for 27–28 March, 24–25 July and 7–8 December 2015 (Figure 3(a)). The correlation coefficients are greater than 0.97, and it demonstrates the calibrated parameters that can be used for estimating regional groundwater flow and the exchange fluxes between the river and the groundwater. The thickness of the sand gravel layer is 20–25 m, and the model parameters, θr, θs, α, n, Ks and l, are 0.065, 0.36, 7.5, 1.89, 20.3 and 0.5, respectively. The thickness of the karst aquifer in the Zhangxia Formation is 140–180 m, and the model parameters are 0.089, 0.25, 0.98, 1.23, 7.293 and 0.5, respectively. The thickness of loam is 20–50 m, and the model parameters are 0.078, 0.43, 3.6, 1.56, 0.2469 and 0.5.
After determining the model parameters, the observed groundwater tables and water discharges during another period of 17 November 2017 to 10 December 2017 were applied for model validation.
(1) Validation with observed groundwater level
Water from the South–North Water Transfer Project has been released into the Yufuhe River since 17 November 2017. The well water table of the S1 monitoring well, which changed considerably from 17 November to 20 December in the initial water release period, was selected as the verification data of the model to reflect the model accuracy well. The correlation coefficient between the measured and simulated water level rise values is 0.97, thereby indicating that the model can simulate the Yufuhe River recharge process well (Figure 3(b)).
(2) Validation with the observed water discharge amount
The water discharge amount reaches the maximum value after 24 days of water release, and the amount of water release is 4.76 million m3. In the modelling, the groundwater level reaches the Quaternary bottom plate when the pressure head of the HYDRUS-3D observation point is 0. The simulated time is 25 days, and the simulated amount of discharge water is 4.93 million m3 (Figure 3(c)). The results show that the model and calibrated parameters are reliable.
The sensitivity analysis is to quantify the uncertainty in the calibrated model caused by uncertainty in the estimates of model parameters and boundary conditions. After sensitivity analysis of the hydraulic conductivity of the karst aquifer (K) and riverbed (C), Manning roughness coefficient (n), and water discharge amount (to change the boundary conditions), the change of boundary condition is the most sensitive to the results. K and n have a similar influence on the effective recharge rate, but their influence on the results are less than that of C.
In accordance with the hydrogeological characteristics of the study area, including the reach lengths and widths, the channel cross-section geometry, energy loss coefficients, and stream junction information, the HEC-RAS model is set up. The streambed is divided into 11 sections (Figure 1) at representative locations throughout the stream reach and at locations where changes in slope, shape, discharge or roughness occur. The critical water depth is selected under the downstream boundary condition, and an implicit finite difference method is used to discretise the calculation. The water stages under the 250,000 m3/d of the water discharge period are simulated, and the results are calibrated by measured data. When the Manning's roughness coefficient of the main channel and the terraces of the river is 0.027 and 0.035, the correlation coefficient between the measured value and the simulated value is greater than 0.95.
Simulation of the effective recharge rate
In the output table of the HEC-RAS results, the loss of surface water flow is 51,840 m3/d. The river stage is transmitted to HYDRUS-3D to simulate the effective recharge rate. When the groundwater depths are 40, 30 and 25 m, the water release times are 26, 20 and 14 days, the water flows are 6.5, 5.0 and 3.5 million m3, the Quaternary losses are 0.934, 0.564 and 0.391 million m3, the surface losses are 1.348, 1.037 and 0.726 million m3, the effective recharge amounts are 4.218, 2.617 and 2.383 million m3, and the effective recharge rates are 64.9%, 65.2% and 68.1%, respectively.
Through the results of the sensitivity analysis, the boundary conditions have the greatest influence on the effective recharge rate. Under the conditions of groundwater depths of 40, 30, and 25 m, water discharge amounts of 350,000 (scenario 1), 300,000 (scenario 2), 200,000 (scenario 3) and 150,000 m3/d (scenario 4) were simulated to analyse the effect of different boundary conditions on the results, and the difference values from the 250,000 m3/d results were calculated respectively. Under scenario 1, the difference values are −13.61%, −11.75% and −14.32% with groundwater depths of 40, 30 and 25 m; under scenario 2, the difference values are −8.48%, −5.41%, −8.31%; under scenario 3, the difference values are 2.30%, 3.08%, 0.59%; and under scenario 4, the difference values are −9.67%, −3.47%, −6.66%. The results show that when the water discharge amount is greater than 250,000 or less than 150,000 m3/d, the effective recharge rate decreases rapidly. The effective recharge rate increases first and then decreases when the water discharge amount is between 200,000 and 250,000 m3/d, which indicates that it is an optimal discharge scheme in terms of the amount of water.
Monitoring data analysis results
The monitoring data analysis method is used with measured data, including surface water and groundwater. During this period, the amount of surface water is 52,037 m3/d obtained by monitoring data and field measurements in the downstream Cuima Bridge section. Furthermore, after analysing the groundwater monitoring data from Quaternary wells Q1 and Q2, the water level, electrical conductivity and temperature suddenly changed when the release water reached monitoring wells given the different quality between the water from the South–North Water Transfer Project and the local groundwater. The time difference between the sudden change points of the two wells is the water movement time in the Quaternary, about 18 days. The average velocity of Quaternary water is 308.9 m/d (the length of stream/the movement time). The area of a mound in the Quaternary is approximately 576 m2 through numerical simulation. The porosity of the Quaternary is 0.18, and the amount of Quaternary water discharge is approximately 32,026.8 m3/d.
On the basis of the statistics of the monitoring wells, the total water discharge during the period is 4.77 million m3. The loss of the surface water is 936,700 m3, and the loss of the Quaternary pore water is 576,500 m3. The sum of the ineffective recharged water is 1,153,200 m3, which accounts for 31.7% of the total water discharge. Therefore, 68.3% of the water infiltrates into the Zhangxia Formation aquifer.
The comparison of the numerical and observed results indicates that the effective recharge rate obtained through the numerical simulation is lower than the observation method. This finding is due to the saturated part below the water surface of the downstream section of the Zhangxia Formation being set as the impervious boundary in establishing the model. In the actual conditions, the karst water of the Zhangxia Formation can flow to the Ordovician aquifer through the Chaomidian fracture which can promote the infiltration of the karst aquifer of the Zhangxia Formation and recharge of the spring. The results of the effective recharge rates of the two methods show that most of the recharged water in the study area can be recharged effectively into the Zhangxia karst aquifer. Thus, the Yufuhe River multi-source water recharging project is necessary. In addition, a high-efficiency and low-cost recharge scheme under different conditions is obtained by analysing the results of the two methods.
The HEC-RAS and HYDRUS-3D semi-coupled model is utilised to simulate the recharge process in various scenarios. When the groundwater depths are 40, 30 and 25 m, the reasonable water discharge periods are 26, 20 and 14 days, and the effective recharge rates are 64.9%, 65.2% and 68.1%, correspondingly. According to the analysis of the conductivity, temperature and water level data of the monitoring well feedback from 17 November 2017 to 10 December 2017, the effective recharge rate is 68.3%.
The results show that the effective recharge rate obtained through the numerical simulation and observation methods is more than 60%, thereby confirming that most of the water can be effectively recharged into the underlying Zhangxia Formation karst aquifer during the recharge process of the Yufuhe River. Through the comparison of effective recharge rates under different water discharge amounts, the effective recharge rate increases first and then decreases when the optimal water discharge amount is between 200,000 and 250,000 m3/d in terms of the amount of water. The results can be used by decision-makers in understanding the amount of recharge water loss and developing a reasonable water release plan. In the next step, a discussion about the transient flow conditions (unsteady and non-uniform flow) and their impacts on the derived results would be further studied and a water quality and quantity dual control model can be established to maintain MAR sustainability.
This study was supported by the Shandong Provincial Key Research and Development Project (Grant No. 2017GSF17121) and the Danish Development Agency (DANIDA) in coordination with the DANIDA Fellowship Centre (Grant No. 17-M08-GEU). The authors would like to acknowledge the editor and anonymous reviewers for their valuable comments, which have greatly improved this paper.