Analytical Solutions of Vertical Airﬂow in an Unconﬁned Aquifer with Rising or Falling Water Table

: The rise and fall of the groundwater level can drive air ﬂow in the vadose zone. In turn, the air ﬂow can interact with the water ﬂow. When the unconﬁned aquifer is covered by a low-permeability media, the coupling of the water ﬂow and the air ﬂow is more obvious. In this study, a conceptual model is established for coupling of air ﬂow and water ﬂow in the vadose zone in response to rapid ﬂuctuations of the water table. Water injection and drainage experiments are conducted in a double-layered sand column with a thick layer (80.5 cm) of coarse sand and a thin layer of ﬁne sand as a low-permeability cap. Different cap thickness (2 cm, 5 cm, 7.5 cm) and different thickness of the vadose zone (30 cm, 40 cm) are set for the experiments. Negative pressure (NP)/positive pressure (PP) of the air in the vadose zone is observed in the drainage/injection experiments, with higher pressure in experiments of thicker cap layer. In each experiment, NP or PP increases rapidly to reach a maximum in the early stage, and gradually becomes zero in the late stage. Analytical solutions on three subdivided stages indicate the permeability and thickness of the cap layer, as well as permeability and porosity of the media in the vadose zone are the key controlling factors on the process of coupling of air ﬂow and water ﬂow. The solutions also reveal the formation mechanism of air pressure in the vadose zone with a low-permeability cap. This study has both theoretical signiﬁcance and engineering applications.


Introduction
The classical theory of groundwater dynamics focuses on saturated flow in aquifers. For groundwater flow in unconfined aquifers, it is generally assumed that the pore air pressure in the unsaturated zone is the same as atmospheric pressure and has nothing to do with groundwater flow in the saturated zone [1][2][3][4][5][6][7]. This assumption works well in many cases but would cause observable errors when the impact of air flow is significant. Vachaud et al. [8] observed that the pressure in the stratified unsaturated zone is significantly different from external atmospheric pressure under a heavy rain. Bouwer and Rice [9] concluded that delayed release of pore water from a pumped, unconfined aquifer can be caused by restricted air movement in the vadose zone. In particular, a low-permeability cap overlying an unconfined aquifer could significantly enhance the difference between pore air pressure and atmospheric pressure driven by water table fluctuation in the aquifer [10,11]. The unsteady-state airflow in a pumped aquifer with falling water table could cause feedback on the change of water table because airflow and groundwater flow are coupled [12,13]. The change of air pressure in the vadose zone driven by falling or rising water table has also been directly observed in experiments with sands [14][15][16][17][18]. Airflow induced by water table fluctuations has important applications in coastal environmental remediation and ecological systems, as well as in unsaturated zone gas transport and natural ventilation of mines and tunnels [10][11][12][13][14][15][16].
In most existing studies of airflow in unconfined aquifers, the coupling between airflow and groundwater flow was considered as a process of multiphase flow, i.e., air-water two-phase flow in unsaturated porous media [10,12,13,15]. To analyze this type of multiphase flow, complex numerical models were built up with a number of parameters [19][20][21][22][23][24][25][26]. The uncertainty of parameter values and scenario-limited results from the models cause difficulty in finding generic behaviors of airflow in an unconfined aquifer. The governing equations of two-phase flow are strongly nonlinear, so analytical solutions incorporating fully the effects of gravity and capillarity in transient multiphase flow through heterogeneous porous media are not tractable [22]. A few of the analytical or semi-analytical models were developed for simplified problems [27][28][29][30][31][32][33][34][35]. Li and Jiao [11] derived an analytical solution for airflow in the unsaturated zone caused by fluctuating water table. Kuang et al. [14] derived nonlinear ordinary equations for changes in air pressure and water table height in a sand column experiment. These studies reveal that simplified models are effective, and suggest that vertical air flow, rather than horizontal air flow, is the primary factor to cause air pressure difference.
This paper focuses on the influence of the vertical water flow on the air pressure difference in the vadose zone, and proposes a simplified nonlinear partial differential equation system. Drainage and injection experiments are conducted in a double-layered sand column with changing thickness of the low-permeability fine sand cap. Analytical solutions are obtained to reveal the general characteristics of the air flow in the case of the water table fluctuations.

Conceptual Model and Assumptions
A conceptual model is established for coupling of air flow and water flow in the vadose zone in response to rapid fluctuation of the water table ( Figure 1). Assume that there are three aquifers of different permeability, with the lowermost aquifer III at a constant water head of z 0 . The high-permeability aquifer III, the medium-permeability aquifer II and the lower part of high-permeability aquifer I are saturated with water. The uppermost aquifer I is capped by a low-permeability fine-sand layer. The vadose zone exists at the upper part of aquifer I and the fine-sand cap layer.  The rise or fall of the water table will drive the air flow in the vadose zone, which  in turn interacts with the water table movement. This phenomenon will be more obvious when the low-permeability cap exists, so that the air cannot flow freely in the vadose zone. As a result, the air pressure in the vadose zone cannot be instantly balanced with the external atmospheric pressure, resulting in a positive or negative pressure relative to atmospheric pressure.

Nonlinear Control Equation of Airflow
The governing equation for vertical one-dimensional airflow is given by generalized Darcy's law for multiphase case where q a is the volumetric flux of air (L −1 T −1 ), k is the intrinsic permeability of the porous medium (L 2 ), k ra is relative permeability of porous medium to air (L 2 ), µ a is dynamics viscosity of air (ML −1 T −1 ), p a is the pressure of air (ML −1 T −2 ), ρ a is the density of air (ML −3 ). Compared to the force caused by air pressure gradients, the force due to gravity is much smaller so it can be neglected under most conditions [36]. Then Equation (1) is linearized to q ma is the mass flux of air (ML −2 T −1 ), which is defined as h a is the air gauge pressure relative to atmospheric pressure in the vadose zone expressed as a height of water column of a reference water density, which is defined as where ρ w is the density of liquid water (ML −3 ), ρ w gh a is relative air pressure in the vadose zone under the fine sand layer, D is the thickness of fine sand layer (L). When the relative pressure of air in the vadose zone is positive, the air in the sand column will be driven upward and the air mass of the vadose zone will be reduced. Conversely, if the vadose zone has a negative pressure, it will attract the air in the atmosphere to flow into the sand column, so that the air mass of the vadose zone will be increased. Assume that the air in the fine sand layer is evenly distributed, and the pore space is not filled with capillary water in the vadose zone. The total air mass in the vadose zone can be expressed as where L is the thickness of the phreatic aquifer (L), h is the thickness of the saturation zone (L), S is the area of the water table (L 2 ), φ is the effective porosity of the coarse sand approximated as a constant. The air density varies with pressure. Considering air as an ideal gas, under the isothermal conditions, according to the Boyle's law ρ = ρ a p/p a , p a is the local atmospheric pressure, and ρ a is the air density under the local atmospheric pressure. The continuity equation of air flow is expressed as According to Equation (3) the air density-pressure relationship, Equation (6) can eventually be rewritten as Equation (7) is the ordinary differential equation indicating air pressure in the vadose zone h a (t). K a reflects the breathability of the fine sand layer.

Feedback Equation for Groundwater Flow
Ignoring the compressibility of porous water and sand, the rate of discharge of water per unit cross-sectional area at any given time is given by Darcy's Law: where q is the rate of discharge per unit cross-sectional area (L 3 T −1 ), K 1 , K 2 is the saturated hydraulic conductivity of unconfined aquifer I and II respectively (L −1 T −1 ), h is the saturated water level of unconfined aquifer at any time. The pressure in the vadose zone is not equal to the atmospheric pressure; the head is not the height of the water table but (h + h a ); h u is the head of the upper plate of the aquifer II; z 0 is the constant head. The rate of recharge or discharge per unit cross-sectional area is equal to the injecting or falling rate of surface of saturation multiplied by effective porosity of the coarse sand [2]: where φ is the effective porosity of the coarse sand, which is approximated as a constant.
In order to solve h, Equating (9) and (10) leads to Equation (11) is the ordinary differential equation indicating the variation of the surface of saturation in the unconfined aquifer.

Experimental Studies
To investigate coupling of air flow and water flow in the vadose zone, we conduct drainage and injection experiments in a double-layered sand column. The sand column in the experiment can be considered as a simplified example of the conceptual model. In the simplified model, only phreatic aquifer is present and the thickness of aquifer II in the conceptual model equals to zero (B = 0). The hydraulic conductivity (K s ) of the phreatic aquifers in the sand column is the same as the hydraulic conductivity of aquifer I (K 1 = K s ) in the conceptual model ( Figure 1). The aquifer III in the conceptual model represents the constant-head reservoir in the experiment, which is a boundary condition with a constant water head. Below we present experimental setup and results that show air pressure in the vadose zone induced by water table fluctuation.

Experimental Setup
The experimental setup (Figures 2 and 3) consists of an organic glass pipe of 104 cm in length, with coarse sand to a height of 80.5 cm inside [37,38]. At the bottom of the sand column, 5 cm of gravel is buried to support and stabilize the water flow. A fine sand layer of 2 cm/5 cm/7.5 cm thickness is placed on the surface of the coarse sand. The saturated permeability of the coarse sand and the fine sand was 70.58 m/d and 4.84 m/d respectively. The porosity of the coarse sand and the fine sand was 0.42 and 0.37 respectively. A U tube is placed 10 cm below the coarse sand surface to measure the pressure difference between the vadose zone and atmosphere. A vent is set under the fine sand layer to keep the pressure in the vadose zone in equilibrium with the atmosphere in the preparation phase of the experiment. During the experiment of water injection or drainage, the vent is closed so that the air in the vadose zone can only pass through the fine sand layer to exchange with the air in the atmosphere.  A constant-head reservoir is used to generate the driving force for the water injection and drainage process. The bottom of the reservoir is connected to the bottom of the sand column through a rubber tube. During the water injection experiment (Figure 2), water is injected into the constant-head reservoir through a high-position glass tube, and the water that overflows from the constant-head reservoir flows into another lower-position glass tube. The change of water level can be converted to the volume of water flowing into the sand column through calculation on water level change of the upper and lower glass tubes. During the drainage experiment (Figure 3), the water flowing through the sand column enters the constant-head reservoir and overflows into a scaled glass tube. The water volume that sand column discharges can be measured by the glass tube. The height of the water level of the constant-head reservoir is z 0 , and the height of the initial water level of the sand column is h 0 . The relative height between the two is reversed in the water injection experiment (z 0 > h 0 ) and the drainage experiment (z 0 < h 0 ). In order to capture the rapid changes of water level in the glass tube, a pressure sensor was put into the glass tube to record the water level automatically.
Before the drainage experiment, the coarse sand and the fine sand were saturated by water for 1 h. Then the water level in the sand column was placed at 20 cm below the surface of the coarse sand, in order to form an initial vadose zone. During the drainage experiment, the constant-head reservoir was lowered to 50 cm (60 cm) below the surface of the coarse sand, so that the head difference was 30 cm (40 cm). The drainage valve was opened at the designated time, and the water levels in the water pipe and the U tube were observed by the pressure sensor and the video recorder. The experiment was continued until the water level in the water pipe was static.
Before the injection experiment, the water level in the sand column was placed at 50 cm (or 60 cm) below the surface of the coarse sand for 1 h. During the injection experiment, the constant-head reservoir was lifted to 20 cm below the surface of the coarse sand, so that the head difference was 30 cm (or 40 cm). Water was supplemented to the constant-head reservoir through the upper glass tube, and the water levels in the water pipe and the U tube were observed by the pressure sensor and the video recorder. The experiment was continued until the water level in the water pipe was static.
Two groups of water injection experiments and three groups of water drainage experiments were performed, corresponding to 2 cm/5 cm fine sand cap in water injection experiments, and 2 cm/5 cm/7.5 cm fine sand cap in water drainage experiments (Table 1) [37,38].

Type Curves of Change in the Air Pressure Beneath the Low-Permeability Cap
The saturated water level cannot be observed by the naked eye, but it can be calculated by the observational result of the cumulative water injection or drainage. Error of the observational data mainly comes from fluctuations in atmospheric pressure and the offset of the sensor. As the amplitude of atmospheric pressure variation is less than 1 cm H 2 O of water column during the experiment, the result of air pressure observation in the vadose zone is reliable.

Results of the Drainage Experiments
During the drainage experiment, significant negative pressure (NP) relative to atmospheric pressure in the vadose zone was observed [37,38]. NP increased rapidly and then gradually returned to zero with time, showing an inverted single-peak curve (Figure 4a). The maximum NP occurred at~10-30 s. The value and the time of maximum NP increase with the increase of fine sand thickness. During the period (10 s to 1000 s), the drainage rate is significantly smaller than in a homogeneous sand column (Figure 4b). The greater the fine sand thickness is, the smaller the drainage rate during this period. NP was not observed in the vadose zone for the homogenous coarse sand column experiment, and was not plotted in Figure 4a. This indicates the coarse sand could not produce the NP in the vadose zone due to the high permeability. Increase of the cumulative outflow in the homogenous sand column was more rapid than those in experiments with a fine sand cap (Figure 4b). This indicates that the fine sand cap promoted to form NP in the vadose zone. The water level of the homogenous sand column fell significantly faster than those in experiments with a fine sand cap, resulting in more than 10 cm difference during the period from 100 s~300 s (Figure 4c). When NP disappeared in the late stage, the water levels in the sand column and the constant-head reservoir gradually became the same. Overall, the experimental results show that significant NP in the vadose zone delays the drop of the water level in the sand column.
From the two sets of experiments with different initial water head differences (h 0 − z 0 ) (Figure 4), it is observed that greater value of (h 0 − z 0 ) corresponds to greater NP and slower rate of water level drop, when thickness of fine sand layer is the same. This is because a certain volume of air was expanded in a short time in the vadose zone, but the atmospheric air could not enter in time via the fine sand layer. The greater amount of expanded air was, the greater the peak NP appeared.

Results of the Injection Experiments
The experimental results show that the positive pressure (PP) in the vadose zone increased rapidly at the beginning and decreased slowly after reaching the maximum [37,38]. The maximum value and its appearance time are related to the thickness of the lowpermeability cap. When the thickness of the fine sand was 2 cm and h 0 − z 0 was 30 cm, the maximum pressure was 10.4 cm of water column at 20 s. When the thickness of the fine sand was 5 cm and h 0 − z 0 was 30 cm, the peak air pressure reached about 17.6 cm of water column at~30−35 s. This indicates that the thicker the low-permeability cap layer was, the higher the peak pressure in the vadose zone reached, and the later the peak value appeared. Corresponding to changes of the air pressure, the water level of the sand column increased gradually from the initial period, increased rapidly in the middle period, and reached stability in the late period.
For the two sets of injection experiments with different initial water level difference (h 0 − z 0 ) but with the same thickness of the fine sand layer (Figure 5), the greater the value of (h 0 − z 0 ) was, the greater the peak positive pressure was. Similar to the drainage experiments, this is because a certain volume of air was condensed in a short time in the vadose zone, but it could not escape to the atmosphere in time via the fine sand layer. The greater the amount of condensed air was, the greater the peak PP appeared.

Limitation Remarks of the Experiments
It is necessary to illustrate that the dynamic model of the water injection and drainage process proposed in this paper does not consider the impact of sand column deformation on the results. During sand injection or drainage, the sand body and Plexiglas tubes undergo a certain change in the stress state, resulting in deformation. When water is drained, the water level of the sand column drops, and there is a vacuum in the vadose zone. At this time, the pressure of the porous fluid decreases, and the effective stress of the sand body increases. The sand column will shrink and deform, thereby strengthening the release of water and weakening the effect of vacuum block. Correspondingly, when the water level rises by the water injection in the sand column, the porous fluid pressure becomes larger and the sand column expands and deforms. However, the shrinking and expanding are too slight to be detected during the experiment, and do not influence the overall experimental results. If the plexiglass sand column is replaced with stainless steel, the disturbance of sand column deformation can be weakened to some extent. In addition, due to the limitations of the experimental conditions in the laboratory, the diameter and height of the sand column used in the experiment are very limited, and the scale is far smaller than the actual distribution characteristics of the formation of the field. When the thickness of the low-permeability cap is large or the distribution is not uniform in actual conditions, the changes of the air pressure may be much more complicated than experiments in the laboratory.

Analytical Solutions
The air-water two-phase equations proposed in this study is basically consistent with the simplified model established by Kuang [15]. This set of equations is nonlinear and it is difficult to obtain a closed analytical solution. The numerical solution can only analyze the overall characteristics, but it is difficult to analyze the relationships between the effect of the air pressure in the vadose zone and the porous media parameters comprehensively. This paper analyzes the dynamic characteristics of the air pressure effect based on the approximate solution of the equation under certain conditions. According to experimental results from this work and previous studies [15,16], three stages are subdivided according to variations of air pressure and surface of saturation (Figure 4c,f). These three stages are resolved respectively.

Early Stage Analytical Solution
In the early stage of the experiment, the drop of the water level was relatively small, h ≈ h 0 . The set of dynamic equations in this case can be simplified as The relative air pressure in the experiments is generally less than 50 cm [15,16], so h a /h a0 < 5%, which is approximated as zero. Then Equation (7) can be rewritten as (14) From this, the approximate solution of the air pressure in the vadose zone of the early stage can be written as C s , C a , K a are constants. C s can be directly calculated based on the basic parameters of the experimental materials. From Equation (15) it can be seen that the early curve develops towards the maximum air pressure C s /C a , and the maximum pressure increases non-linearly with increasing thickness of fine sand layer thickness D.

Middle Stage Analytical Solution
In the middle stage of the experiments, the rate of drainage or injection is becoming steady, and dh/dt can be regarded as a constant, dh/dt = C. The airflow equation in this case can be simplified as The airflow equation in this case can be simplified as From this, the approximate solution of the air pressure in the vadose zone of the middle stage can be written as where h δ is an integral constant. At this stage, the maximum air pressure in the vadose zone is positively correlated with the rate of change of the water level in the unconfined aquifer.

Late Stage Analytical Solution
In the late stage of the experiments, the drainage process becomes very slow, the position of the water surface is almost unchanged, there are approximately dh/dt = 0 and h ≈ z − h a , the airflow equation in this case can be simplified as Thus, the characteristics of pressure change in the later stage of experiments are obtained.
where h c is an integral constant, which has the same unit as the pressure. Since the air pressure also decays to a small value in the late-stage experiments, h a /(L − z 0 ) can be approximated as zero, so the logarithmic approximation of the pressure is linear with time.
From the ln-h a -t curve, K a /D can be derived.

Estimating Parameters with Semi-Analytical Solutions
The conceptual model and the analytical solutions can be widely used under natural conditions. The drainage/injection experiment can be considered as a simplified example of the conceptual model (B = 0, K 1 = K 2 = K s ). The dynamic equations and the dynamic analysis features in early and late stages of the experiments show that the pressure in the vadose zone is controlled by the two parameters of K s / φ and K a /D. K s / φ is included in the C s parameter of Equation (16), which reflects the influence of the high-permeability layer (coarse sand) on the drainage and injection process. K a /D is included in the parameter C a of Equation (16), which reflects the influence of the low-permeability layer (fine sand) on the drainage or injection process. K s / φ = 0.299 cm/s can be calculated according to the properties of the experimental materials. C s = 9.89 cm/s can be calculated from Equation (17). In contrast, K a /D is an uncertain parameter.
The analytical solutions and the measured data in the first 30 s of the drainage experiment were fitted. K a /D was obtained from the results of the drainage experiments at 2 cm, 5 cm, and 7.5 cm thickness of the fine sand layer. The values of K a /D are 2.3 cm/s, 1.1 cm/s and 0.6 cm/s respectively (Figure 6a), which are generally inversely proportional to the thickness of fine sand layer. The corresponding K a is 4.5~5.5 cm 2 /s, a parameter reflecting the air permeability of the fine sand layer in the early stage of experiment. C a can be derived from Equation (16) as 0.28~0.36 s −1 . The maximum pressure is mainly determined by C s /C a , and the maximum value tends to increase nonlinearly with the increase of the fine sand layer thickness D. When D→∞, the peak pressure will reach the limit value (h 0 − z 0 ), which is the difference in the original water head. According to Equation (25), the observed data for the drainage time in the range of 1000~3000 s and the analytical resolutions were fitted. The data for the final period (>3000 s) was not used due to the little change and poor accuracy of the U-tube. Based on the experimental results with fine sand thickness of 2 cm, 5 cm, and 7.5 cm, the fitting results for K a /D are 12.9 cm/s, 9.5 cm/s, and 8.6 m/s respectively (Figure 6b), approximately inversely proportional to the thickness of the fine sand. The corresponding K a is 26~65 cm 2 /s, a parameter that reflects breathability of the fine sand layer in the late stage of drainage.
Based on the analysis above, it can be inferred that the K a calculated from the early drainage data is one order of magnitude lower than that derived from the late drainage data, which indicates that the relative air permeability value of fine sand layer is not constant during the experiment. The reason for this is that in the fine sand there was residual moisture which weakened the breathability of the fine sand layer. In the late stage, the long-term drainage decreased the moisture content of the fine sand and increased the breathability.

The Generic Behaviors of Airflow Revealed by the Analytical Solutions
Overall, it is indicated from the analytical solutions that in the early stage, the maximum air pressure develops towards C s /C a , and the maximum pressure non-linearly increases with the fine sand layer thickness D. In the middle stage, the maximum air pressure in the vadose zone is positively correlated with the rate of the water level change in the unconfined aquifer. In the late stage, the logarithmic approximation of the pressure is linear with time, and K a /D can be derived from the ln h a versus t curve.

Conclusions
Changes in the water table affect air flow in the vadose zone. When the unconfined aquifer is covered by the low permeability layer, the air will not flow freely and will form a considerable positive pressure or negative pressure in the vadose zone relative to atmospheric pressure. In this study, a simplified nonlinear model describing the experimental phenomena is established. Water drainage and injection experiments conducted in a double-layer sand column with a low-permeability fine sand cap layer are analyzed. In the drainage experiment, as the water level of the sand column decreases, the negative pressure in the vadose decreases rapidly firstly, reaches a maximum value and then increases slowly to zero. In the water injection experiment, as the water level of the sand column rises, the positive pressure in the vadose zone increases rapidly and then decreases slowly to zero. The thicker the fine sand layer is, the greater the NP or PP in the drainage or water injection process is, and the longer time is needed to reach the maximum value. It is verified experimentally that the maximum vacuum degree increases nonlinearly with the increase of fine sand cap thickness. Analytical solutions on three subdivided stages are provided. Only three parameters are used to explain the changes of the air pressure in the vadose zone and the saturated water level in the dynamic model. The air pressure in the vadose zone mainly depends on the saturated permeability of the coarse sand, the breathability of the fine sand, and its thickness. The air pressure in the vadose zone in the early and late drainage stages can be described by an approximate analytical formula.