Evaluation model coupling exploitable groundwater resources and land subsidence control in regional loose sediments

The loose sediments in the Yangtze River Delta, the North China Plain, the plain of Northern Jiangsu and other districts in China are of great thickness, complex in structure and abundant in groundwater. Groundwater overexploitation easily results in geological disasters of land subsidence. Aiming at the issues, assessment models coupling exploitable groundwater resources and land subsidence control in regional loose sediments were brought up in this paper. The two models were: (1) a three dimensional groundwater seepage model with land subsidence based on the one dimensional Terzaghi consolidation theory; (2) a three dimensional full coupling model on groundwater seepage and land subsidence based on the Biot consolidation theory to simulate and calculate. It can be used to simulate and calculate the problems in real situations. Thus, the groundwater seepage and land subsidence were coupled together in the model to evaluate the amount of exploitable groundwater under the specific requirements of land subsidence control. The full coupling model, which considers the non-linear characteristics of soil mass and the dynamic changes of soil permeability with stress state based on the Biot consolidation theory, is more coincident with the variation characteristics of the hydraulic and mechanical properties of soil mass during the pumping process, making the evaluation results more scientific and reasonable.


Introduction
The Quaternary and Tertiary loose sediments in the Yangtze River Delta, the North China Plain, the Northern Jiangsu Plain and some other districts in China are of great thickness, which preserves abundant groundwater resources. The sediments consist of Holocene pore-phreatic aquifers (group), late Pleistocene No. Ⅰ confined aquifers (group), mid-Pleistocene No. Ⅱ confined aquifers (group) and early Pleistocene No. Ⅲ In this study, based on groundwater seepage, the one dimensional Terzaghi consolidation theory and the Biot consolidation theory, a three dimensional groundwater seepage model coupling land subsidence and a three dimensional full coupling model on groundwater seepage and land subsidence were brought up. The groundwater seepage and land subsidence were coupled together in the model to evaluate the exploitable groundwater resources under the requirements of land subsidence control. The full coupling model considered the non-linear characteristics of the soil mass and the dynamic changes  are the cosine of the exonormal of the discharge flow rate boundary's and the coordinate axis's included angle;   is the saturation deficiency (an increase of free surface) or specific yield (free surface is down), which indicates the water yield that unit section area of aquifer absorbs or discharges when free surface changes a unit height ;  w q is the recharge of atmospheric precipitation on the unit area of the free surface (m/d);  (Luo Zu-jiang and Zhang Yue-ping 2007). The calculation model of the deformation extent of the aquifer caused by the decline of the groundwater level is as follows: The elastic deformation of the confined aquifer is The inelastic deformation of the confined aquifer is The elastic deformation of the phreatic aquifer is The inelastic deformation of the phreatic aquifer is There into:  In the subsidence model, the aquifers and the aquitards inside were treated as a whole. The whole skeleton of the aquifer was assumed to be compressed and potentially resulted in land subsidence. The parameters fv S and fv S were equivalent across the whole scale of the aquifer. The total amount of the land subsidence is calculated only after the amount of compression in every aquifer is obtained.
There were two steps for the coupling of the flow model, land subsidence model and the water release quantity produced by the deformation of soil layers. First, the water level should be worked out using the flow model and take it as a boundary condition of the subsidence model. The deformation amount of soil layers and the resulting water release quantity can be calculated and substituted into the flow model. Finally, compute them with repeated iterations to finish the coupling of water flow and land subsidence.
By coupling the hydraulic heads in the two models together, the coupling simulation model of groundwater exploitation and land subsidence was established.

The Three Dimensional Full Coupling Model on Groundwater Seepage and Land Subsidence
The change of the groundwater seepage field caused by groundwater exploitation is actually a process of interaction and mutual effects between the seepage field and the stress field. Based on the Biot consolidation theory, the three dimensional full coupling model on groundwater seepage and land subsidence has considered the non-linear characteristics of soil mass and the dynamic changes of soil permeability with stress state, making the results closer to reality and more precise.
Where:  G is the shear modulus;   is the Poisson's ratio;  x w , y w and z w are the displacement components of directions x , y and z;  u is the pore water pressure;   is the density of earth's gravity;  w  is the density of water's gravity. Combined with fixed solution conditions, such as the initial conditions and boundary conditions, the equation above can be solved by the numerical method.

Fixed Solution Conditions.
 The Initial Conditions A. The Initial Condition of Ground Stress Using the dead-weight stress of soil to estimate its initial stress, Where:  C and  are intensity-indexes, and k, n, Rf, G, F and D are test parameters achieved by the conventional triaxial compression test of soil.
 Dynamic Model of Permeability k The fluid-solid coupling is actually the stress redistribution of the soil skeleton caused by dissipation of pore water stress, which leads to the consolidation deformation on the macro level and changes of porosity factor n in soil. Therefore, the permeability of the soil is changed and the seepage is also consequently affected. At present, many models have been proposed in papers about fluid-solid coupling ignore the interrelationship between porosity and permeability and the changes of them. Under the premise of the assumptions of the Biot consolidation theory, the dynamic expression of porosity φ and permeability k was built as follows: In which:  0 k is the initial intrinsic permeability;  V  is the volumetric strain. The model above couples groundwater seepage with soil deformation in its mechanism. The model involves a large amount of parameters, and the solving of the model is consequently complex. Therefore, at present it is still rare to use the model in regional evaluation of groundwater resources.

Instance One
3.1.1. General Situation, Generalization and Discretization of the Study Zone. Jiangyin City is located in the Yangtze River Delta economic zone, on the south of the Yangtze River, facing Jingjiang across the river, on the north of the Shanghai-Nanjing expressway, on the west of Changshu and Zhangjiagang and on the east of Changzhou and Wujin. It is at the geometric centre of the "Golden Triangle" of Suzhou, Wuxi and Changzhou. The pore water-bearing system in Quaternary shallow loose rock within 60 meters contains pore-phreatic aquifer, micro-confined pore aquifer, clay aquitard and No.Ⅰconfined aquifer from top to bottom. The internal structure of the system was generalized as inhomogeneous and anisotropic in three dimensional space. The Yangtze River was treated as a constant water head boundary, and the rest of the boundaries around the system were treated as general water head boundaries. The top of the system was treated as a free surface boundary (Bridget R Scanlon et al. 2002), and the bottom of the system was treated as an impervious boundary. The groundwater exploitation took towns as units and was treated as big wells.
The three dimensional model above coupling groundwater seepage and land subsidence based on the one dimensional Terzaghi consolidation theory was solved by the finite difference method

Identification and Validation of Parameters.
Using the model above, the test lasted for 5.25 days according to the data of the observation wells. The whole time, the segment was discretized into 15 stress periods and each stress period was divided into 10 steps. The initial flow field of each aquifer was measured in the test, and the initial displacement was set as 0. The pumpage of each pumping well was obtained by actual statistics. The hydraulic conductivity on the boundaries and the initial parameters of the subareas in each aquifer were all based on previous data, field tests and indoor experiments, and the water level of observation wells were fitted with pumping test data. Then, the inverse calculation method was used to compute the hydro-geological parameters. The parameters involved in coupling analysis included: initial horizontal permeability coefficient K0x and K0y, initial vertical permeability coefficient K0z, specific yield μ, internal friction angle φ, cohesion force c, initial elastic modulus E0, initial poisson ratio v0 and unit weight γ, among which the soil mechanics parameters were obtained by indoor experiments. Figure 3 and table 1 illustrate the relevant parameter division in the 4th layer of the model and the parameter values in each subarea. Figure 4 and figure 5 illustrate the fitting precision of the model of groundwater level and land subsidence, from which we can see that the fitting precision of the groundwater level is relatively high and the land subsidence is fitted to the subsidence trend. Thus, the model can be used to simulate the change of water level and land subsidence caused by groundwater exploitation in the study zone.   K0y The 4th Confined Aquifer  Ⅰconfined aquifer, which would not be lower than half the depth of the roof of the aquifer in 10 years, and the quantity of the land subsidence, which would not exceed 50mm in 10 years, as constrained conditions. After the identification of the model, the groundwater recoverable resources of towns in Jiangyin and the mining well number were obtained by running the model, and the results are listed in table 2. Figure 6 and figure 7 are the predictive groundwater level contour maps of the micro-confined aquifer and the 1st confined aquifer at the end of April 2017, and figure 8 is the predictive land subsidence contour map at the end of April 2017. The result shows that the shallow groundwater exploitation in Jiangyin is of great potential as long as the exploitation layout of the shallow groundwater is planned rationally. The total amount of exploitable resources in the micro-confined pore aquifer is 1492.16×10 4 m 3 /a, the total amount of exploitable resource in No.Ⅰconfined aquifer is 2052.47×10 4 m 3 /a and the maximum quantity of the land subsidence will not exceed 50mm in 10 years.

General Situation, Generalization and Discretization of the Study Zone.
Wujiang City lies in the alluvial plain of the Yangtze River Delta, and it has deposited a huge thickness of Quaternary loose sediment. The pore water-bearing system in Quaternary shallow loose rock buried within 60 meters contains phreatic aquifer, clay aquitard and micro-confined aquifer from top to bottom. The internal structure of the system can be generalized as inhomogeneous and anisotropic in three dimensional space. The top of the system was treated as a free surface boundary, and the bottom of the system was treated as an impervious boundary. Considering the exploitation of shallow groundwater in the surrounding cities, the study zone surroundings were generalized as impervious boundaries. The initial water level and initial displacement were given by statistical material. The groundwater exploitation took rural towns as units and was treated as big wells.

Identification and Validation of Parameters.
The model above was simulated using the Groundwater Seepage and Land Subsidence Software (GWS) that was developed by the research group. According to the data of the observation wells, the test lasted for 1.25 days. The whole time segment was discretized into 10 stress periods, and each stress period was divided into 10 steps. The initial flow field of each aquifer was measured in the test, and the initial displacement was set as 0.
The pumpage of each pumping well was obtained by actual statistics. The hydraulic conductivity on the boundaries and the initial parameters of the subareas in each aquifer were all based on previous data, field tests and indoor experiment, and the water level of the observation wells were fitted with pumping test data. Then, the inverse calculation method was used to compute the hydro-geological parameters. The parameters involved in the coupling analysis included: initial horizontal permeability coefficient K0x and K0y, initial vertical permeability coefficient K0z, specific yield μ, internal friction angle φ, cohesion force c, initial elastic modulus E0, initial poisson ratio v0 and unit weight γ, among which the soil mechanics parameters were obtained by indoor experiments. Figure 10 and table 3 illustrated the relevant parameter division in the 1st aquifer of the model and the parameter values in each subarea. Figure 11 illustrated the fitting precision of the model of the groundwater level and land subsidence, from which we can see the fitting precision of the groundwater level is relatively high. Thus, the model can be used to simulate the change of the water level and land subsidence caused by groundwater exploitation in the study zone.  Figure 10. Parameters of the first confined aquifer.  Figure 11. Fittings of groundwater level and subsidence with data from the simulation and observation.

Prediction of Groundwater Level and Land Subsidence.
The model regarded the end of April 2007 as the starting point of the calculations, and the buried depth of the groundwater level in the micro-confined aquifer, which would not be lower than half the depth of the roof of the aquifer in 10 years, and the quantity of the land subsidence, which would not exceed 50mm in 10 years, as constrained conditions. After the identification of the model, the groundwater recoverable resources of towns in Wujiang and the mining well number were obtained by running the model, which are listed in table 4. Figure 12 and figure 13 are the predictive groundwater level contour maps of the phreatic aquifer and the micro-confined aquifer at the end of April 2017, and figure 14 is the predictive land subsidence contour map at the end of April 2017. The result shows that the shallow groundwater exploitation in Wujiang has great potential as long as the exploitation layout of the shallow groundwater is planned rationally. The total amount of exploitable resource in the micro-confined pore aquifer is 821.50×10 4 m 3 /a, the total amount of exploitable resource in the phreatic aquifer is 382.09 ×10 4 m 3 /a and the maximum quantity of land subsidence will not exceed 50mm in 10 years.

Conclusions
The evaluation model coupling exploitable groundwater resources and land subsidence control can calculate the exploitable groundwater resources based on the land subsidence control in regional loose sediments. The model couples the exploitable groundwater resources with the land subsidence control, thus unifying and integrating the groundwater exploitation and land subsidence control effectively. Based on the Biot consolidation theory, the full coupling model considered the non-linear characteristics of the soil and the dynamic changes of soil permeability with the stress state. It depicted the variation characteristics of the hydraulic and mechanical properties of the soil mass during the pumping process more correctly than the three dimensional model coupling groundwater seepage and land subsidence based on the one dimensional Terzaghi consolidation theory, and the result of the evaluation is more scientific and reasonable. However, the model is not yet mature enough to be widely used in the evaluation of the exploitable groundwater resources and land subsidence control in regional loose sediments at present because it has too many parameters and the equation is difficult to solve. With the continuing development of calculation technology and the accumulation of the data of the groundwater regime and the stress-strain of the soil with the water, it is inevitable that the theory and method will be widely used in the future.