Abstract
Two analytical solutions using segregation variable method to calculate the hydraulic head under steady and unsteady flow conditions based on Tóth’s classical model were developed. The impacts of anisotropy ratio, hydraulic conductivity (K), and specific yield () on the flow patterns were analyzed. It was found that the area of the equal velocity region increases and the penetrating depth of the flow system decreases at steady state with anisotropy ratio increases, which is defined as . In addition, stagnant zones can be found in the flow field where the streamlines have opposite directions. These stagnant zones move toward the surface as the horizontal hydraulic conductivity increases. The results of the study on transient flow indicate that a relative increase in hydraulic conductivity produces a positive impact on hydraulic head and a relative enhancement in specific yield produces a negative effect on hydraulic head at early times.
1. Introduction
Tóth [1, 2] analytically calculated patterns of gravity-driven groundwater flow systems for small drainage basins using the hydraulic potential theory of Hubbert [3]. It was found that only one large regional groundwater flow system could be developed for a homogeneous and isotropic basin with a linearly sloping water table. However, local, intermediate, and regional groundwater flow systems could develop with a water table attributable to the topographic relief of a sinusoidally undulating water table superimposed on a linearly sloping straight line. The nested groundwater flow systems are essential for understanding a variety of geologic processes, such as fluid flow in sedimentary basins [4], surface water-groundwater interaction [5], regional heat and solute transport [6], and the role of groundwater in the study of ore deposition and petroleum accumulation [7, 8]. Tóth [9] presented an extensive overview of the history, principles, methods, applications, and natural effects of gravity-driven groundwater flow system.
Furthermore, the physical mechanisms of groundwater flow system have been explained by Tóth [4] and Engelen and Jones [10]. Numerical simulations, including the concepts and methods [11–18] and afterword physical simulations [19, 20], have also been developed. The theory of groundwater flow systems has been applied to study many practical problems in science and engineering [21–24].
A progression of groundwater flow system investigations has been conducted based on the assumption that the medium is isotropic and homogeneous [2, 25, 26] or heterogeneous [14, 27, 28]. For example, Tóth [2] investigated the impact of different depth models on the groundwater flow systems; and Jiang et al. [14] investigated the impact of hydraulic conductivity on groundwater flow systems by using the COMSOL software.
Based on a critical review of the literature, one can see that most of the previous studies focused on the steady-state flow of groundwater. However, the groundwater flow is usually transient in the fields. Models of steady-state groundwater flow can result in significant limitations in applications. Besides, the prevailing analytical solutions of unsteady groundwater flow are essentially targeted on determining well flow issues. There are few specific reports about analytical solutions for unsteady nested groundwater flow systems [29]. In order to have a better understanding about the evolution process of groundwater flow systems from unsteady to steady under natural and mining conditions, it is necessary to carry out the theoretical investigation on the unsteady groundwater flow systems. This study focuses on the analytical methodology of the two-dimensional (2D) unsteady groundwater flow systems in Tóth’s classical model, which will provide new insights into understanding the groundwater flow systems.
2. Hydraulic Head Solution under Steady-State Condition
In most of the previous studies, the media were assumed to be homogeneous and isotropic [25, 26, 30], and the heterogeneity has been considered in few studies [27, 28]. However, the porous medium is always anisotropic. The horizontal hydraulic conductivity is usually greater than the vertical one. In this paper, an analytical solution of the steady-state groundwater flow system in a small basin considering the anisotropy of hydraulic conductivity based on Tóth’s classical model [2] (Figure 1) was derived.
The governing equation of the potential function for a 2D steady-state flow field considering the effect of anisotropy can be described as follows:
The boundary conditions can be written asin which (m) is the hydraulic head (); (m/d) and (m/d) are the horizontal and vertical hydraulic conductivity, respectively; (m) is the length of the basin; (m) is the depth of the basin, and the datum level is chosen at the bottom of the basin; represents the average slope of the valley flank; and represent the amplitude and frequency of the sinusoid, respectively; , , and represents the wavelength.
The following dimensionless variables are defined: , , , , , and . So the problem can be transformed to the following equations:
The mathematical model can be solved by using the segregation variable method. The details can be found in Appendix A. The solution of the hydraulic head can be expressed aswhere ,
If we set , , , , , , , and in (4), (4) will be changed to
Equation (6) is the same as the solution of Tóth [2]. One can use the new solution equation (4) to investigate the impact of anisotropy on groundwater flow systems. The new solution can also be used to study the unsteady flow. A MATLAB program [31] was developed to calculate (4). The parameters were given as m, m, , and . Then , . According to Darcy’s law,
The dimensionless hydraulic head, , and the velocity distribution were calculated with those parameters and different values of anisotropy, , as shown in Figure 2.
(a)
(b)
(c)
The area of equal velocity region decreases significantly with the anisotropy ratio (the contour interval is 0.003). The penetrating depth of groundwater flow systems decreases with the horizontal hydraulic conductivity, or the horizontal flow velocity increases with the horizontal hydraulic conductivity. Additionally, stagnant zones can be found from the simulating results. They are located at the places that have opposite direction for the streamlines (stagnant zone-1) and at both sides at the bottom in the field (stagnant zone-2). The stagnant zone-1 moves toward ground surface as increases.
3. Hydraulic Head Solution under Unsteady-State Condition
3.1. Mathematical Model and Its Analytical Solutions
According to Wang and Wan [29], a 2D unsteady flow may be considered as a steady flow plus a relative head. The initial condition is ( is a constant). Assuming , where is the hydraulic head at the steady flow (4), is relative head. Thus, the hydraulic head at unsteady flow can be obtained by solving and then plus . The mathematical model can be described as follows:
The initial condition and boundary conditions can be written aswhere (m−1) is specific yield and (d) is the time.
Similarly, some additional dimensionless variables can be defined: , , , and . Thus, the problem can be transformed to the subsequent equations:
The potential function can be obtained using segregation variable technique (the detailed derivation can be found in Appendix B):where
Equation (15) is the same as (4) when the time is large enough (, in Figure 4).
3.2. Results and Discussions
A MATLAB program was developed to calculate (15). The parameters are given as m, m, , , m, m−1, and m/d. So , , , , and . When time is large enough, the flow approaches quasi-steady state, which means that the hydraulic head does not change with time (Figure 3(c)).
(a)
(b)
(c)
The hydraulic heads at the right side are significantly higher than those at the left side if a sinusoidally undulating water table is set at (Figure 3). In addition, the distribution of is irregular at early times (Figure 3(a)). When time is large enough, the ultimate distribution of unsteady flow is similar to that of the steady state (Figures 3(c) and 2(a)). Figure 4 shows the -time behavior for three different positions A, B, and C (Figure 1). Each point desires a precise time to achieve a steady state indicating the flow from unsteady to steady, and the flow at the right side approaches steady state faster than that of the left side.
3.3. Sensitivity Analysis of the Parameters
Sensitivity analysis is a tool to analyze the impact of the input parameters on the results of a model. It is essential to analyze the sensitivity of parameters in order to improve the accuracy of groundwater flow models and to reduce errors induced from the uncertainty of hydrogeological parameters [32]. Studying the response of simulation results to the variations of important parameters can not only facilitate improving the reliability and accuracy of the models but also play a vital role in hydrogeological survey.
In unsteady groundwater flow systems, hydraulic head and specific yield can have great impacts on the evolution process of unsteady flow. The normalized sensitivity method [33, 34] was used to analyze the sensitivity of the results to and . The normalized sensitivity of ) and to the relative modification of a given parameter can be expressed as follows:where is of point changes with time before the parameters change and is of point changes with time after the parameters change.
The values of and were increased by 10%, 50%, and 100% to investigate their impacts on the hydraulic head at the middle of the domain (, Figure 1) (Figure 5). Figure 5 shows that a relative increase in produces a positive effect on hydraulic head and a relative increase in produces a negative effect on hydraulic head at early times. Hydraulic head is sensitive to and at early times but insensitive to them at late times.
4. Summary and Conclusions
On the basis of Tóth’s classical model, analytical solution of hydraulic head containing hydraulic conductivity under steady and unsteady flow conditions is obtained. From this study, the following conclusions can be drawn:(1)For the steady flow, the area of equal velocity region becomes much smaller with a larger anisotropy ratio, and the penetrating depth of groundwater flow systems becomes smaller. Stagnant zones locate at the places which have opposite directions for the streamlines and both sides at the bottom of the field. Finally, when becomes larger, stagnant zones at the places with opposite streamlines become much closer to the surface.(2)For the transient flow, when time is large enough, the ultimate distribution of unsteady flow is consistent with the steady groundwater flow system model. The closer it is to the right side, the faster it is for the flow to approach steady state.(3)For the transient flow, a relative increase (e.g., 10%, 50%, and 100%) in produces a positive impact on hydraulic head, and a relative increase in produces a negative impact on hydraulic head at early times.
Appendices
A. Analytical Solution of Potential Function for Steady-State Flow
Analytical solution of hydraulic head for (3) is presented as follows.
Let
One can obtain two ordinary differential equations:
According to the left and right boundary conditions,
Considering the bottom boundary condition,
Discussions. Consider the following(I)When , it does not meet the actual situation.(II)When , general solution for (A.2) is . Using the condition (A.4), . General solution for (A.3) is . Using the condition (A.5), . Let , so .(III)When , let ( is nonzero real number).
General solution for (A.2) is
According to boundary conditions, one can get that
General solution for (A.3) is
Using the condition (A.5),
Therefore,
Combined with the upper boundary condition, one can obtain(A.12) is the Fourier series; then
Then one can obtain , .
B. Analytical Solution of Potential Function for Transient Flow
Analytical solution of hydraulic head for (10) is presented as follows.
Let
One can obtain three ordinary differential equations:
According to the boundary conditions, one can get
Discussion of . Consider the following(I)When , it does not meet the actual situation.(II)When , general solution for (B.2) is . Using the condition (B.5), .(III)When , let ( is nonzero real number).
General solution for (B.2) is
Therefore,Summing up the situations of and ,
Discussion of . When , it does not meet the actual situation, so , letting .
General solution for (B.2) is .
Therefore,
General solution for (B.3) isSo,
Derived from the principle of superposition,According to the initial condition,
So, The left of (B.17) and (B.18) is the Fourier series of the right side; then one can get by using the method of integration by parts. Then one can obtain
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This research is sponsored by National Natural Science Foundation of China (Grants nos. 41272258, 41372253, 41521001, and U1403282), National Basic Research Program of 973 Program (Grant no. 2010CB428802), and the Fundamental Research Funds for the Central Universities, China University of Geosciences (Wuhan) (Grant no. CUG140503). The authors also would like to thank the editor and the anonymous reviewer for rendering valuable comments and suggestions for improving this paper.