Geofluids

Volume 2018, Article ID 1054730, 19 pages

https://doi.org/10.1155/2018/1054730

## Quantify the Pore Water Velocity Distribution by a Celerity Function

^{1}Key Laboratory of Meteorological Disaster, Ministry of Education / Joint International Research Laboratory of Climate and Environment Change / Collaborative Innovation Centre on Forecast and Evaluation of Meteorological Disasters / School of Hydrology and Water Resources, Nanjing University of Information Science and Technology, Nanjing, Jiangsu 210044, China^{2}Water Resources Section, Faculty of Civil Engineering and Geosciences, Delft University of Technology, Stevinweg 1, PO Box 5048, 2600 GA Delft, Netherlands^{3}Department of Physical Geography and Geoecology, Faculty of Science, Charles University, Albertov 6, 128 43 Prague 2, Czech Republic^{4}Key Laboratory of Mountain Hazards and Surface Process, Institute of Mountain Hazards and Environment, Chinese Academy of Sciences, Chengdu 610041, China

Correspondence should be addressed to Ye Su; zc.inuc.rutan@us.ey

Received 2 February 2018; Revised 4 April 2018; Accepted 7 May 2018; Published 3 July 2018

Academic Editor: Joaquín Jiménez-Martínez

Copyright © 2018 Wei Shao et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

Celerity and velocity of subsurface flow in soil porous medium are intimately linked with tracer transport, while only few current studies focused on their relations. This study conducted theoretical analyses based on a pore bundle model, under which celerity in unsaturated flow is equivalent to maximum velocity. Furthermore, under 3 models Brooks-Corey model and unmodified and modified Mualem-van Genuchten models, we used a celerity function to derive breakthrough curves aiming to quantify the advective tracer transport for 5 typical soil textures. The results showed that the celerity under near-saturated conditions can be 5 up to 100 times larger than the saturated hydraulic conductivity, and a small volumetric fraction (<15%) of pores contributed more than half of the specific discharge. First arrival time and extensive tailing of the breakthrough curves were controlled by maximum velocity and velocity distribution, accordingly. As the kinematic ratio in the Brooks-Corey model remained constantly for each specific soil, we used it to quantify the ratio of maximum tracer velocity over average tracer velocity. We also found that a bimodal soil hydraulic function (for a dual-permeability model) may result in similar soil hydraulic conductivity functions for different parameter sets, but their celerities are different. The results showed that the proposed celerity function can assist in investigating subsurface flow and tracer transport, and the kinematic ratio could be used to predict the first arrival time of a conservative tracer.

#### 1. Introduction

Subsurface flow is conducted in the soil porous medium, in which the pore structure (i.e., pore connectivity, tortuosity, and pore size distribution) controls the permeabilities and pore water velocities in numerous flow paths [1–3]. In hydrological studies, it is impractical to deterministically represent the pore-scale fluid dynamics and transport processes in soil, for the variability of pore water velocity is difficult to observe or predict [4–6]. For describing the subsurface flow processes, the continuum approach was often used to conceptualize the discrete materials and water particles as one single continuous mass. The soil hydraulic properties (e.g., the lumped hydraulic conductivity and moisture capacity) and state variables (e.g., volumetric water content and capillary pressure) are often quantified in a representative elementary volume (REV) [7, 8]. The continuum subsurface flow equation (i.e., the well-known Darcy-Richards equation) calculates specific discharge and average pore water velocity for transport phenomenon. Specifically, the specific discharge quantifies the volumetric flow [9, 10], and the average pore water velocity can be coupled with an advection-diffusion equation to simulate solute/thermal transport [10–12]. Applications of the Darcy-Richards equation relay on the indicators of velocities and specific discharge.

For describing the water fluid dynamics, the terminologies of velocity (i.e., pore water velocity) and celerity have been commonly used in various hydrological systems, such as river channels [13–15], estuaries [16, 17], and soil porous medium [8, 18, 19]. Theoretically, velocity expresses mass flux, while celerity represents the speed of pressure propagation through a flow domain [20]. In surface water, for instance, the pressure wave is a flood wave propagating through river channels. Similarly, the celerity in subsurface flow describes a perturbation-induced pressure wave that is affected by precipitation, evaporation, or fluid injection and extraction [21, 22]. In saturate soil, the difference between velocity and celerity can be illustrated by a virtual experimental [20], and the authors pointed out: “*in a cylinder full of sand and saturated with water, changing the flow rate or head at the input boundary will immediately cause a change in flow at the output boundary. While the water flow velocity through the sand is slow, the celerity in this case is (theoretically) instant, hence the immediate response. At larger scales, this case is analogous to a confined aquifer with incompressible water and rock. Allowing for the compressibility will slow the celerity a little, but the velocities of flow will still be much less than the celerities.*”

In confined saturated flow, the nearly-instant celerity is caused by the low compressibility of water and soil porous medium. While in unsaturated soil, the pressure propagation was achieved by the variation of capillary pressure that is driven by different mechanisms. For example, a response of the pore water pressure in unsaturated soil can be caused by either preferential flow [23] or pressure wave through the entrapped air [24]. However, how to distinguish the initial cause-effect mechanism of pressure response is still a challenge [22, 25]. Many studies used the kinematic wave equation to simulate unsaturated flow [8, 18, 26, 27], in which the celerity was defined as a first-order derivative of hydraulic conductivity (i.e., equal to specific discharge) with respect to the water content [22, 28, 29]. These studies illustrated the celerity as “mobility of water in soil” [29, 30], and the mechanisms of pressure propagation in unsaturated zone was, however, not discussed.

On the other hand, the pore water velocity in numerous flow paths is dictated by the pore structure such as pore connectivity, tortuosity, and pore size distribution [1–3]. Especially, the soil hydraulic functions can be derived as an integration of pore water velocities in all flow paths if considers a unit gradient condition [10, 12, 31]. Based on a pore bundle model, the hydraulic properties of each flow path can be equivalent to that of a cylinder tube with a small circular cross section and a sufficient length [32–34], then the experimentally measured water retention curve is used to infer pore size distribution and pore water velocity distribution [35–38]. Under such condition, each tube conducts the viscous flow; therefore, the relations of the equivalent tube radius to capillarity pressure and pore hydraulic conductivity can be, respectively, described by a capillary rise equation (i.e., Young-Laplace equation) and a pipe/cylinder liquid flow equation (i.e., the Hagen-Poiseuille law). The hydraulic conductivity functions integrated from pore velocity distribution can indirectly reflect the pore-scale hydraulic properties [39, 40]. As the first-order derivative of hydraulic conductivity can be defined as a celerity function for unsaturated soil [22, 24], the celerity can be related to the pore water velocity distribution.

The complex solute transport behaviours are caused by the pore water velocity distribution in a subsurface flow system [28, 41–43], which have been captured either by an experimentally obtained breakthrough curve [5, 44–46] or by natural and artificial tracers in a hydrological system [28, 47, 48]. The breakthrough curve expresses the time series of the tracer concentration of drainage flow, which was often characterized by a coexisting early initial breakthrough and extensive tailing [44]. The tracer experiment in a hillslope or a catchment assists in inferring the travel time distribution [49, 50]. However, many tracer studies found a paradox behaviour; that is, during a high-intensity rainstorm a large amount of “old” water resident in soil or groundwater is expelled into the stream, while the labelled “new” water of the rainfall on the hillslope immediately starts to appear in the stream after a fast infiltration through a fast flow path [51–54]. Those complex tracer transport behaviours address the importance of quantifying the pore water velocity distribution for achieving the more detailed physical interpretations of unsaturated flow and transport [28, 55].

The celerity and velocity are intimately linked with tracer transport, while most of the current studies focused on one or two aspects. For example, mathematical derivations for the celerity and kinematic ratio have been provided by Rasmussen et al. [22], while they did not find the relation between tracer transport and pressure propagation. Wang et al. [24] and Mohammadi et al. [56] found the derivation of the tracer breakthrough curve through the soil hydraulic conductivity functions, and their method have been validated through their experiments. However, the studies from Wang et al. [24] and Mohammadi et al. [56] were based on the soil hydraulic functions of either the Brooks-Corey model or the van Genuchten model for a single-permeability model, and more soil hydraulic functions still need further studies. There is no report existing that explored the relations among celerity (pressure propagation), velocity (water transport), and tracer transport (pore velocity distribution) yet.

The objective of this paper was to (i) illustrate the concept of celerity in unsaturated flow, (ii) propose a new soil hydraulic function (i.e., celerity function) to express the pore water velocity distribution, and (iii) quantify the advective breakthrough curves using the proposed celerity function. Specifically, in Section 2.1, the definitions of velocity and celerity were given to illustrate the mechanisms of pressure propagation in unsaturated soils. Based on a pore bundle model, Section 2.2 mathematically proved the equivalence between celerity in unsaturated soil and maximum pore water velocity among all flow paths that conduct water. Furthermore, a theoretical analysis was conducted to illustrate the difference between celerity and velocity of wetting front under a unit-hydraulic-gradient condition, showing that the celerity function derived from a perturbation analysis can reveal the pore water velocity distribution (see Section 2.3). The celerity functions of the Brooks-Corey model and both unmodified and modified Mualem-van Genuchten models are provided in Section 2.4. The impact of velocity distribution on tracer transport can be illustrated by mathematically deriving a breakthrough curve through a celerity function (Section 2.5). Functions for the celerity in dual-permeability models are derived in Section 2.6. Section 3 analysed the hydraulic characteristic and breakthrough curves of typical soil textures by using the proposed celerity function. Hereto, the pore water velocity distribution and breakthrough curve of typical soils were analysed using the existing parameter sets of the Mualem-van Genuchten model and the Brooks-Corey model. Additionally, the equifinal parameter sets of a bimodal soil hydraulic function were distinguished by the celerity function. The follow-up discussion focused on the mechanism of pressure propagation and tracer transport in a natural subsurface flow system.

#### 2. Theory

##### 2.1. Definitions

In subsurface, the vertical component (LT^{−1}) of a specific discharge vector (i.e., the volume flux of water per unit cross-sectional area) is formulated using Darcy’s law as
where (LT^{−1}) is the hydraulic conductivity, (L) is the pressure head, and (L) is the vertical coordinate (positive downward).

The average pore water velocity (LT^{−1}) is
where (−) is the volumetric water content and (−) is the residual water content.

The continuity equation for one-dimensional vertical flow is

The derivative of specific discharge with respect to is written as

Substitution of (4) for in (3) results in an advection equation
where (LT^{−1}) is the celerity [22]

Theoretically, the advection equation represents the advection of the soil-water content with speed . When is a constant value, it would mean that an arbitrarily shaped pulse of moves with constant speed without changing its shape. The celerity is, however, not a constant but a function of water content. Hence, the celerity is an approximate advection speed of a small pulse of .

A ratio between celerity and average velocity is called the kinematic ratio (−), as defined by Rasmussen et al. [22]:

##### 2.2. Celerity and Maximum Pore Velocity

The pore sizes and their intrinsic permeabilities vary throughout the unsaturated zone. Based on a pore bundle model, many previous studies considered soil as a bundle of nonintersecting, parallel, and cylindrical tubes with varying radii to quantify the water and solute transport [32–34, 57]. When water enters soil flowing through a pore of certain size, each pore is either entirely filled with water or entirely empty. Considering variably saturated conditions, pores are filled from the smallest tube group () to the largest (). Under unsaturated conditions, of the tube groups are filled with water. The pore water velocity (LT^{−1}) in one specific equivalent (saturated) cylinder tube can be determined by its specific discharge
where (LT^{−1}) is a coefficient relating the head gradient and the average velocity in pore .

The pores with the same size and the same pore water velocity can be classified as one group, based on which the water content (i.e., all the volumetric fractions of the pore space filled with water) is discretised in fragments of equal water content and ordered from small pores to larger pores. Assuming a unique, nonhysteretic relationship between water content and capillary pressure, water can flow through the water-filled saturated smaller pores with higher capillary pressure, while the larger pores will remain empty and inactivated for conducting unsaturated flow [58].

By discretizing the flow domain into numerous flow tubes that are vertically homogeneous, the hydraulic characteristics of each tube can be formulated as functions of the volumetric water content . Considering the pore size distribution, the pore water flow velocity follows a distribution ordering from smallest pore group (with velocity ) to largest pore group (with velocity ). At a certain water content value, (LT^{−1}) may denote the maximum pore water flow velocity taking place in the water-filled pore with the largest equivalent radius—then the specific discharge (LT^{−1}) and the hydraulic conductivity (LT^{−1}) can be obtained through a summation:
where is the largest tube group filled with water. Substitution of (8) for in (9) gives
where the hydraulic conductivity is defined as

The unit hydraulic gradient condition considers a gravity-driven vertical infiltrating flow (i.e., when and ), and the hydraulic conductivity function therefore can be obtained by summing up the pore water velocity from a limited number of pore groups as shown in Figure 1.