2 Basic Relationships and the Soil Phases Under Equilibrium Conditions

Chapter 2: Basic Relationships and the Soil Phases Under Equilibrium Conditions


2.1      Where we are going

2.1.1    Soil as a Three-Phase System

2.1.2    Volume and Mass Relationships

2.1.3    Liquid Content

2.2       Characterization of the Solid Phase

2.2.1    Texture and systems of classification

2.2.2    Sieving

2.2.3    Sedimentation

2.2.4    Surface Area of Soil Particles

2.2.5    Soil Structure and Bulk Density

Aggregation, bio-pores, packing geometries

2.3       Characterization of the Liquid Phase

2.3.1    Soil Water Content Measurement Methods

Gravimetric, dielectric, neutron probe, cosmic ray backscatter, remote sensing via radar and gravitation (SMAP, SMOS, GRACE); direct, indirect, and surrogate methods, errors.

2.3.2    Applications of Soil Water Content Information – The Water Balance across scales of soil profile, field, catchment, continent, and globe and seconds to decades.

2.3.3    Field Capacity, Wilting Point and Plant-Available Soil Water

2.4       The Energy State of Soil Water

2.4.1    Total Soil Water Potential

2.4.2    Gravitational Potential

2.4.3    Matric Potential Water Properties Affecting the Matric Potential

2.4.4    Pressure Potential

2.4.5    Solute Potential

2.4.6    Water Potentials Under Equilibrium Conditions (Hanks, 1992 – pp. 31-37)          Constructing Potential Diagrams under Equilibrium Conditions          Water Potentials Under Nonequilibrium Conditions

2.5       Bringing gas, solids and liquids together

2.5.1    Capillarity and Surface Tension What is surface tension Deriving the LaPlace Equation Temperature dependence of surface tension Surface absorption (Tuller and Or Equation) and the myth of dry soil Contact angle

Idealized, impacts of surface texture and defects, hysteresis

2.5.2    Soil Water Content-Energy Relationships Soil Water Characteristic  (SWC) Measurement of SWC Curves Fitting Parametric SWC Expressions to Measured Data Hysteresis in the Soil Water Characteristic Relation          The independent domain model for hysteresis          Analytical expressions for hysteresis

2.5.3    Measurement of Soil Water Potential Components Tensiometer for Soil Matric Potential Measurement Heat Dissipation in Rigid Porous Matrix for Measuring Matric Potential The Psychrometric equation and Psychrometric methods for Measuring Water Potential Piezometer for Measuring Hydrostatic Pressure Potential Measurement of Soil Water Osmotic Potential

2.6       SUMMARY

2.7       Review Questions

2.8       REFERENCES

2.9       PROBLEMS

2.10     Appendices:

2.10.1  The thermodynamics of liquid and solid (Hamaker) interfaces

2.10.2  Advanced methods in soil moisture measurement (GPR, ERT, EMI, NMR, COSMOS

2.10.3  The independent domain model for hysteresis

2.1 Where we are going

We are concerned with describing the retention and movement of fluids in natural porous media. This objective fits neatly into a single sentence, however the challenge of attaining this goal remains to be met to the degree required to solve many contemporary problems in environmental sciences and engineering. Nonetheless, there are many powerful tools at our disposal. In this chapter we will cover the basic physical properties upon which much of the remainder of the text rests. If all goes well you will get a “feel” for both the physical processes which control flow through porous media and the conceptual and mathematical approaches we employ to characterize the processes. Having a well balanced understanding from both perspectives is essential; you must be able to guide your application of mathematical solutions based on a clear and accurate intuitive concept of what you are trying to describe.How do you describe a chunk of moist porous medium? What are some basic features which would allow you to abstract that chunk of material in terms of a set of mathematical parameters which will assist in describing: the way in which it might dry; how rain falling on it would soak in (and how much would run off); and what would happen to it if you spilled some gasoline on it? We require a system of measurable quantities for a porous medium that may be used in addressing the broadest possible set of questions about how fluids will interact with it under a variety of conditions. Porous media is characterized by a set of parameters selected to provide as concisely as possible the greatest insight into the response of the media to a range of physical processes. In this section we identify the fundamental parameters employed in describing both the media and the fluid of interest, and introduce the fundamental equations of motion for fluids in porous media.

2.1.1 Soil as a Three-Phase System

Before digging into the details of hydraulic properties and flow processes, we need to provide some background, definitions and quantitative measures. We divide natural porous media into three constituent phases: solid, liquid, and gaseous, with each of these phases consisting of an admixture of compounds.

The gaseous phase is generally dominated by the major constituents of the atmosphere; N2, O2, CO2, H2O vapor etc. Clearly, as microbes and plant roots respire, a soil’s CO2 and CH4 content may become significantly elevated, while in industrially contaminated sites there may be significant levels of organic vapors. By and large we keep track of the state of the gaseous phase in the soil just as we do in the atmosphere; we consider the partial pressures of each of the constituents, and keep track of the partitioning of each compound. Primary parameters of interest, then, include the gas diffusivity, and the partial pressure versus temperature relationship, and the Henry’s law coefficient (Henry’s law states that at a given temperature, for low concentrations, there is a linear relationship between the partial pressure of the gas and the concentration of the material in the liquid phase). It is noteworthy that transport of gaseous phases are dominated by molecular diffusion which generates about 2 cm/day of spreading, while molecular diffusion in the liquid phase is 100 times less effective in contaminant mass movement. The movement of gaseous phase constituents is also a function of the changes in atmospheric pressure, which generates a pumping effect, as noted by Buckingham (1904) and Weeks (1979). Of course, gaseous phase movement can also be driven by liquid phase movement, as is the case of a saturated wetting front pushing gas down through a profile during flood irrigation. Gas phase transport has received renewed interest with the advent of soil vapor extraction (SVE) as a remediation technique for volatile organics. In this procedure, gas is pumped out of the vadose zone through a well, taking the volatile fraction of contamination with it.

Of the possible liquid phase constituents, water is by far of greatest interest. In the vadose zone it is valid to assume liquids are incompressible, since by definition liquids will only experience pressures which are negative, and generally with a magnitude of less than 10 atmospheres. In contrast, below the water table, groundwater commonly experiences extreme pressures (e.g. in cases of confined and very deep aquifers). Otherwise liquid phase properties are most often parameterized in terms of the interaction with solid and gaseous phases, such as contact angle, capillary pressure, volatility, etc., as is discussed later.


What is the percent volumetric change in water in the vadose zone with a change in pressure from zero to –10 bars (10 bars is 106 N/m2)? Recall that the isothermal compressibility of water, bw, is about 4.8 x 10-10 m2/N. The percentage change in volume of water is 100 x bw x the change in pressure = 0.048%. Since the material is mostly gas filled when the pressure is this low, this very small change in volume will be undetectable.

The solid phase may be as complicated or simple as you like. Typically composed of largely inert minerals (rock, gravel, sand, silt and clay), solid constituents may include roots, carbon, and other highly reactive substances. Note that in this discussion we have slighted the semi-solid or soft tissue “goo” such as microbes, plants and worms that are prevalent in the vadose zone. These components are very important to understanding the chemical behavior of the unsaturated zone, and are discussed in chapter 4.

To get yourself in the proper frame of mind, pretend you are an ant. What humans consider to be uniform porous medium now looks lumpy and irregular. Although this may seem all too obvious, it is crucial to consider this problem in detail. To do this we will have to introduce the basic parameters used to conceptualize a porous media. First, let’s identify the constituents which we will concern ourselves with: the mass and volume of each phase; gas, liquid and solid. As identified in Figure 2.1a, the volume and mass of phase i are identified as Vi and Mi respectively.

Figure 2.1a. Diagram of constituents of a sample of unsaturated porous media on volumetric and mass basis. Notice that the air and water phases are interchangeable on short time scales (down to seconds), while the other constituents are largely static over timescales of months and years.

2.1.2 Volume and Mass Relationships

We can now define the dry bulk density of the medium, which is the dry mass per unit volume of the medium, is denoted ρvb.

(1)   \begin{align*} \mathrm{\uprho_{vb}} &= \frac{\text{mass in volume V}} {\text{volume V}} \nonumber \\ \nonumber \\ &= \mathrm{\frac{M_{vs}} {V_{t}}} \end{align*}

Here the subscripts “Vt” denotes the total volume over which this density is calculated, “b” denotes bulk and “s” denotes the solid phase (without fluid). Typical units for density are grams per cubic centimeter, kilograms per cubic meter, or, heaven forbid, pounds per cubic foot. The upper limit on ρvb is the density of the pure solid mineral, often called the particle density, which we denote ρs.

(2)   \begin{align*} \mathrm{\uprho_{s}} &= \frac{\text{mass of pure mineral}} {\text{volume of pure mineral}} \nonumber \\ \nonumber \\ &= \mathrm{\frac{M_{s}} {V_{s}}} \end{align*}

Notice that the value of ρs is independent of the value of Vs (since the value of density is taken to be a material property independent of the volume selected), whereas as we shall soon see, ρvb is not quite as cooperative. The solid phase density is often denoted by the specific gravity of the material, denoted by s, which is simply the ratio of the solid density to that of water ρw:

(3)   \begin{equation*} \mathrm{s=\uprho_{s}/\uprho_{w}} \end{equation*}

Typical values of s are 2.65 for quartz (often used as the default value in lieu of further information), 2.54 for feldspar, 2.72 for calcite and 5.0 for pyrite.

Complementing these definitions, is that of porosity, or void fraction, which we denote by nv.

(4)   \begin{align*} \mathrm{n_{v}} &= \frac{\text{sample volume} - \text{volume of pure material}} {\text{sample volume}} & &\text{(a)} \nonumber \\ \nonumber \\ &= \mathrm{\frac{V_{t} - V_{s}}{V_{t}}} & &\text{(b)} \nonumber \\ \nonumber \\ &= \frac{\text{void volume}} {\text{sample volume}} & &\text{(c)} \nonumber \\ \nonumber \\ &= \mathrm{1 - \frac{V_{s}} {V_{t}}} & &\text{(d)} \nonumber \\ \nonumber \\ &= \mathrm{1 - \frac{\uprho_{vb}} {\uprho_{s}}} & &\text{(e)} \end{align*}

nv indicates the fraction of the volume V which is not filled by solids. The bulk density of a medium, ρvb, depends primarily on two characteristics of the medium: its particle size distribution and the medium’s structure. Figure 2.1 illustrates that the size of particles alone does not dictate density: A set of media which are made up of single sizes of particles (and which are packed in the same arrangement) will have identical densities. Hence, the particle size distribution, and to a lesser degree the packing, strongly affect ρvb. The material structure of the medium is also significant in determining ρvb, particularly in the case of clays and very sharply angular aggregates.

So far so good; but what about this annoying V subscript? Can’t we just forget about the sample volume? Not really, and it is worthwhile seeing why. Let’s just run a quick “experiment” to make this clear. Imagine that you are in the middle of a large field with a tool which measures ρvb for any prescribed V. Start with the instrument set with V smaller than any individual grain in the soil (Figure 2.2). You will be sampling the density of either pure solid, pure void, or maybe some arbitrary mixture of the two. Now, holding the instrument steady, enlarge V continuously, all the while recording the values of ρvb. Figure 2.2 shows two possible realizations of such an experiment.

Figure 2.1
Figure 2.1 Dependence of bulk density on particle size distribution : (a) uniform particle size distribution gives uniform density for a given geometry of particle; (b) increasing the range of particle sizes gives rise to greater bulk density.  Note that “well sorted” means that the particles fall within a narrow range of sizes, which is the opposite of “well graded,” which refers to media containing a wide range of particle sizes (or “grades”).
Figure 2.2 Bulk density as a function of control volume size.  The Representative Elementary Volume (REV) for porosity is attained at the scale identified as v3.  Note that the REV is not unique for any given media, but rather is a function of the parameter which is being characterized.  Thus, since permeability depends on the largest connected apertures, which might be very widely spaced, the REV for permeability might be much larger than that found for porosity.

Although the measured density becomes “better behaved” as the control volume enlarges from the grain scale, there is no scale over which density is entirely constant. As the control volume increases, ever larger structural features of the setting are encompassed with differing properties, and thus the measured bulk density changes. The existence of such multi scale heterogeneity is widely accepted, and Dagan (1986) provides a very nice treatise on the topic. Keep this in mind: the value of all physical parameters in the vadose zone are a function of both position and sample volume.

We have two choices then: 1) throw up our arms, and figure that natural systems are hopelessly complex, or, 2) make some reasonable simplifying assumptions. For the sake of progress, we will opt for the latter option (however tempting the former). We will assume that there is some volume, much larger than the grains of the porous medium, yet smaller than the distance between dissimilar regions, which provides a representative sample of our porous medium. This volume is referred to as the Representative Elementary Volume (REV) for the system (Bear, 1972). Each extensive property of the medium is defined at a point by its value calculated taking V to be one REV about that point. We do get around to a more complete analysis of this problem in chapter 6 where we explicitly address scale dependent variability.

With that settled, we can proceed to complete our laundry list of definitions of physical properties of porous media (dropping the v subscripts since this is now implicitly understood):

Total Bulk Density : ρt (solid and liquid phase)

(5)   \begin{equation*} \mathrm{\uprho_{t} = \frac{M_{t}} {V_{t}}} \end{equation*}

where Mt is the total mass (solids and liquids).

Void Ratio: e (ratio of pore volume to solid volume)

(6)   \begin{equation*} \mathrm{e} = \mathrm{\frac{V_{a} + V_{w}} {V_{s}}} \end{equation*}

where Va is the air volume. Note that using our earlier definitions, we may relate the void ratio to the porosity:

(7)   \begin{align*} \mathrm{e} &= \mathrm{\frac{V_{a} + V_{w} + V_{s} - V_{s}} {V_{s}}} & &\text{(a)} \nonumber \\ \nonumber \\ &= \mathrm{\frac{V_{t} - V_{s}} {V_{s}}} & &\text{(b)} \nonumber \\ \nonumber \\ &= \mathrm{\frac{1 - V_{s} / V_{t}} {V_{s} / V_{t}}} & &\text{(c)} \nonumber \\ \nonumber \\ &= \mathrm{\frac{n} {1 - n}} & &\text{(d)} \end{align*}


Table 2.1: Range of values found in bulk density and associated properties illustrated by particular cases. Particle density is taken as 2650 kg m-3.

Description Bulk Density, rb

kg m-3

Porosity, n Void Ratio, e
Surface of wet clay soil 1120 0.58 1.37
Surface of loam soil 1280 0.52 1.07
Spheres (uniform size) open packing 1390 0.48 0.91
Subsoil of sandy texture 1610 0.39 0.65
Compacted sandy loam soil 1900 0.28 0.39
Spheres (uniform size) close packing 1960 0.26 0.35
Sandstone 2120 0.20 0.25

(Source: Marshall and Holmes, 1988; with permission of Cambridge University Press)

2.1.3 Liquid Content

Central to most investigations of vadose zone processes is the quantification of liquid content (typically water, often referred to as the moisture content). This is expressed in either a mass or volumetric basis, as shown in 8 and 9, with the volumetric being more commonly employed, either as a water content or degree of saturation.

a. Mass basis (a.k.a., gravametric)

(8)   \begin{equation*} \mathrm{w = \frac{M_{w}} {M_{s}} = \frac{r_{w} V_{w}} {r_{s} V_{s}} = \frac{r_{w} V_{w}} {r_{b} V_{b}}} \end{equation*}

b. Volume basis. θ is called the volumetric liquid content and S the degree of saturation:

(9)   \begin{equation*} \mathrm{\uptheta = \frac{V_{w}} {V_{t}} = w \frac{\uprho_{b}} {\uprho_{w}}} \end{equation*}

(10)   \begin{equation*} \mathrm{S = \frac{\uptheta - \uptheta_{0}} {\uptheta_{s} - \uptheta_{0}}} \end{equation*}


θo = “residual” or “irreducible” liquid content.

θs = saturated liquid content.

Equivalent Depth of Soil Water (Soil Water Storage): De = θD, where D is the soil depth increment having water content θ. De provides a means for expressing water content in a given soil depth increment as soil water storage (Length), and for relating aboveground water dimensions (Length; as in rainfall, irrigation, or Et) to belowground dimensions (soil water content, L3 water/L3 soil). This is very useful in water balance calculations, for example.

The typical units for moisture content are in mass per mass (for gravimetric water content, e.g., g/g) or volume per volume (for volumetric water content, e.g., cm3/cm3). Watch out though, because you will see water content in other units akin to the depth of storage; how about inches of water per foot of depth? Typically the units will indicate whether the values are on a volume or mass basis.

We need to take some time to iron out a few necessary fussy points regarding θs and θo. Lets start with the easy one: what could be more straightforward than the concept of saturation? Just fill all the pores with water, right? Well, although this is simple in principle, it is rarely achieved in the vadose zone. Due to the complexity of the pore space, it is common for some pores to remain gas filled even if there is ponded water on the surface due. This is due to trapping in dead end pores and pores that become completely surrounded by water filled pores. For this reason, several terms have been proposed to describe the state of natural soil being as wet as it will get under normal conditions. This is referred to as a “satuated”, “satiated”, or “field saturated” condition. This can be used in the definition of degree of saturation if clearly specified. Typically, the field saturated conditions render 90% of the pore space filled with fluid (Bower, 1966).

Now, for this business of residual moisture content. Unless soil is heated to over 700 oC , water held in hydrogen bonds will remain in place (Nagata et al., 1974, Figure 2.XX). This is clearly an extreme which will not be reached in typical natural settings, but the point is that the residual moisture content a particular medium will achieve is a function of the drying process. There are three drying processes which are widely considered in the definition of θo: gravity drainage, plant uptake, and atmospheric surface drying. The residual moisture content achieved under gravity drainage is referred to as “field capacity,” θfc, which is defined as the condition which arises after a soil has been allowed to drain to a stable condition under gravity driven conditions (assuming an infinitely deep water table, for a detailed discussion, see Assouline, S., and D. Or, 2014, The concept of field capacity revisited: Defining intrinsic static and dynamic criteria for soil internal drainage dynamics, Water Resour. Res., 50, 47874802, doi:10.1002/2014WR015475). From an operational perspective, this is often assumed to be the conditions which occurs when the soil water is at a pressure of -1/3 bar, although clearly such a threshold is very much a function of the soil texture, ranging from -1/30 bar for coarse sands to -1 bar for clayey soils (Assouline et al., 2014). With more time, water eventually stops draining when the films of water break between soil particles leaving isolated pockets of water.

Figure 2.XX Experimental mass versus drying temperature for sepiolite clay, illustrating that 3% of the mass is made up of bound water beyond 700℃ (with permission, Nagata et al., 1974). [Nagata, Hroshi, Susumu Shimoda, and Toshio Sudo. “On dehydration of bound water of sepiolite.” Clays & Clay Minerals. 1974.]

The second commonly taken value for residual moisture content is referred to as the “permanent wilting point,” θpwp. You might have guessed that this is the moisture content which is observed when plants can no longer retrieve water from soil (always lower than θfc). Since plants vary greatly in their ability to extract water, this is operationally defined as the moisture content when the water in the medium is a pressure of -15 bar. Since water no longer flows along continuous streams at moisture contents less than θfc this measurement relies on vapor transport of water between pores, and so is necessarily slow in achieving equilibrium. The final method of drying is due to direct evaporation to the atmosphere. This process can bring the upper soil surface to pressures of less than -100 Bar. In summary, so long as you define the conditions over which a flow process will occur, there will be a real and measurable residual moisture content, but if you are unsure of the context where the parameter will be used, it is more reasonable to provide a range of possible values, or provide the extreme values just discussed. In fact, many soil scientists believe that θo should be simply taken as 0, regardless of the problem under study, as this is the only physically defensible value for arbitrary conditions. Table 2.1a and Figure 2.X provide some typical ranges for various parameters.

Figure 2.X Typical values of saturation, field capacity and permanent wilting point for a range of soil textures. The Plant Available Soil Water (PASW) indicates the volume fraction of the soil which might retain water following infiltration and drainage (i.e., field capacity), and be possible for the plant to extract (i.e., up to the PWP).

A primary practical importance of the field capacity and wilting point concepts is for determination of a range of plant-available soil water (PASW; Fig.2.X). Soil water storage available for plant use is generally calculated as being between these two limits (θFCWP), as water contents higher than θFC are generally not sustained for long times. It takes only one to a few days for the attainment of field capacity after wetting if conditions permit. Plant available soil water storage is an important factor in the determination of irrigation amounts for a cropped field or other soil-plant system. For practical purposes irrigation amounts in excess of a soil’s field capacity are lost to deep percolation and thus should be avoided.

A useful ‘rule of thumb’ in absence of other information is to estimate θFC as θS/2, and θWP as θFC/2; i.e., a soil exhibiting this property will have lost one-half it’s saturated water content at field capacity, and another one-half the remaining water by wilting point.

Knowledge of soil water content is also important to geotechnical applications where the attainment of high bulk density by compaction minimizes future settlement and increases the soil strength. For each level of compactive effort there is an “optimal” water content at which maximum compaction is attained (this occurs at an approximate degree of saturation of 0.8). Hence, knowledge of the soil water content can assist a geotechnical engineer in determining how much water should be applied to a layer of earth fill to achieve optimal compaction. Such information is also useful for farmers who want to avoid excessive compaction of their fields by farm implements (see Hillel 1980, Chapter 14 for details).

Table 2.1a Typical values of physical properties for several soils.

Soil Texture † Porosity, n
* Saturated Hydraulic Conductivity, Ks
† Specific Density, s † Field Capacity
†Permanent wilting point
Sandy 0.38 10-3 – 10-1 1.65 0.09 0.04
Sandy Loam 0.43 10-3 – 10-2 1.50 0.14 0.06
Loam 0.47 10-3 – 10-2 1.40 0.22 0.10
Clay Loam 0.49 10-4 – 10-3 1.35 0.27 0.13
Silty Clay 0.51 10-7 – 10-4 1.30 0.31 0.15
Clay 0.53 0 – 10-4 1.25 0.35 0.17

†Adapted from Israelson and Hanson, Irrigation Principles and Practices, Wiley and Sons, 1962.
*May vary beyond the stated range.

Example 1-1): Volume-Mass Relationships

Problem Statement:

The diameter of a cylindrical soil sample is 0.1 m and its height is 0.3 m. (1) What was the gravimetric water content (W) if the wet mass was 4.12 kg and the dry mass was 3.18 kg ? (2) What was the soil bulk density, ρb? (3) What was the volumetric water content, θv ?; and (4) What was the porosity (n), the void ratio (e), and the degree of saturation (S) ?


(1) Gravimetric water content W

(11)   \begin{equation*} W = \frac{m_{w}}{m_{s}} = \frac{4.12 - 3.18}{3.18} = 0.295 \left[\frac{\si{kg}}{\si{kg}} \right] \end{equation*}

The gravimetric water content W is the mass of water relative to the mass of oven-dry soil.

(2) Soil bulk density ρb

(12)   \begin{equation*} \rho_{b} = \frac{m_{s}}{V_{t}} = \frac{3.18}{2.35 \cdot 10^{-3}} = 1353 \left[ \frac{\si{kg}}{\si{\cubic\meter}} \right] \end{equation*}

The dry soil bulk density rb is the ratio of the mass of oven dry soil to its bulk volume. The soil bulk volume is equal to the volume of the cylinder given as:

(13)   \begin{equation*} V = A \cdot h = r^{2} \pi h = 0.05^{2} \cdot \pi \cdot 0.3 = 2.35 \cdot 10^{-3} \left[ \si{\cubic\meter} \right] \end{equation*}

(3) Volumetric water content θv

(14)   \begin{equation*} \theta_{v} = \theta_{m} \frac{\rho_{b}}{\rho_{w}} = 0.295 * \frac{1353}{1000} = 0.399 \cong 0.4 \left[ \frac{\si{\cubic\meter}}{\si{\cubic\meter}} \right] \end{equation*}

The volumetric water content θv is the volume of water relative to the soil’s total volume and may be expressed as gravimetric water content W times the ratio of bulk density ρb to water density ρw = 1000 kg/m3

(4) Porosity n, void ratio e, and degree of saturation S

(15)   \begin{align*} &n = 1 - \frac{\rho_{b}}{\rho_{s}} = 1 - \frac{1353}{2650} = 0.49 \\ &e = \frac{n}{1 - n} = \frac{0.49}{1 - 0.49} = 0.96 \\ &S = \frac{\theta_{v}}{n} = \frac{0.399}{0.489} = 0.81 = 81\% \end{align*}

The porosity n is the ratio of pore volume to the total volume and may be expressed as 1 minus the ratio of bulk density to particle density. The particle density is assumed to be 2650 kg/m3

2.4       The Energy State of Soil Water

Imagine a column of saturated soil that is allowed to drain: water will flow downward draining the upper soil, but why?  What drives this motion? Water content alone is not sufficient to characterize the water status of soils. An additional critical attribute is the energy state of soil water. Soil water contains different forms and quantities of energy, as does any body of matter. Some familiar forms of energy are kineticpotential, thermal, electrostatic, and chemical. The last aspect of characterizing the state of water in porous media requires the introduction of the energetic status of the liquid in the porous medium. While all of these quantities represent energy, we will  discover that most liquid movement in porous media results from gradient in only the mechanical potentials (as opposed to electrostatic, chemical or thermal potential). In principle this is as simple as water flowing downhill in a stream: water is flowing from higher elevation (thus higher gravitational potential) to lower elevation under constant pressure (the gauge pressure of water at the air-water surface is zero), thus moving from areas of high gravitational potential to areas of lower potential with constant pressure potential.  But given the differing roles of the different potentials, we must be deliberate in identifying these terms and their role in soil-water processes.

Let’s start with the potential that you experience every time you move: Kinetic energy (Ek) is acquired by virtue of motion and is proportional to the velocity squared \mathrm{[E<sub>k</sub>=1/2 mv^{2}]}. Illustrating how different the primary processes are in natural porous media compared to daily life, because the movement of water in soils is generally quite slow (e.g., of the order of meters per year), its kinetic energy is almost always negligible. We call these slow-moving flows “creeping flow” to indicate that the forces of acceleration associated with changing kinetic energy are much smaller than the forces generated by viscosity and gradients in pressure.  These flows are associated with Reynolds numbers less than 1.  For every general rule there are special case exceptions, and in soils this can here occur during extreme events: flow along macropores can have velocities with kinetic energy that are appreciable. We will discuss this in detail in the context of dynamics of fluid in porous media, but for the moment we simply posit that potential energy is the primary driver of flow outside of these rather extreme cases.

With respect to potential energy we need to discriminate between terms which are mechanical from those which are electrostatic, chemical or thermal. The point here is that although hot water has more energy than cold water, this energetic difference does not cause the hot water to flow toward the cold water if they are in physical contact (think of a lake in summer with warm water on top and cool water below) . Similarly, a pure water has higher osmotic potential than does a saline water, but when you have salty water below fresh water you will not see the salty water flow in bulk toward the fresh. So how can these chemical or thermal potentials give rise to movement of water?  To translate chemical potential into mechanical potential requires there be an interface which selectively constrains the motion of selected components of the system (i.e., an charge barrier, permeable membrane, or gas-liquid interfaces). Some examples are helpful here. Consider a semi-permeable membranes that allows water to pass but excluding solutes.  This can result in water delivery to the side with higher solute concentration, resulting in bulk water movement.  At a clay surface with electrostatic bound cations population, when exposed to water, these cations will hydrate, pulling in bulk water into the clay.  Systems with such selective regions can develop huge mechanical pressures (e.g., turgor in cells or the pressure drop across reverse-osmosis desalinization systems). Such interfaces occur in the natural world at roots, cell walls, and air-water interfaces (Figure 2.X123). While very important to life, these processes are somewhat esoteric in the case of the movement of bulk water in soils.

Figure 2.X123. The total water potential arises due to multiple mechanisms, including: (a) Pressure potential which can result from surface tension, where curved air-water interfaces give rise to pressure following the La Place equation, or if pure water is held in a rigid vessel with a semi-permeable membrane, the pressure can become very negative if the vessel is in a low-humidity environment; (b) Gouy-layer potential arising from water interactions at mineral surfaces, which generally includes contributions from dipole-dipole attraction of water to the mineral, electrostatic attraction of the positively charged portion of water to the negatively charged mineral surface, electrostatic attraction of the negatively charged portion of the water molecule to the cations bound to the mineral surface, and osmotic potential acting on the water as it dilutes the ions attracted to the minerals; and (c) the osmotic potential that can attract water vapor at saline water surfaces. (a)-(c) taken together are known as the matric potential (to the degree this term has a universally agreed upon definition), and may be expressed in units of pressure or head. When the gravitational potential is added this is known as the total water potential.

To predict fluid movement in the vadose zone we are generally concerned with keeping track of the total mechanical potential (TMP which we will indicate with the letter H) of the fluid phase. The total mechanical potential is made up of the sum of the gravitational potential (due to elevation) and the pressure potential (the force per unit  area exerted by the fluid on its boundary – Figure 2.X124). As we introduced above, fluid motion results from a gradient in the TMP.  So we are concerned with gradients in the TMP.  Still we need to define reference (or datum) in order to provide a frame of reference  from which the TMP might vary.  For pressure we typically define the atmospheric pressure as the datum (by assigning it a value of 0). This is the pressure of the water table.  In selecting a datum for gravity, we choose an elevation which is convenient for each problem. In the lab this may be the elevation of the surface a pool of water, and in the field it is often taken to be the soil surface. We emphasize that the specific choice of datum is immaterial: fluid moves relative to local changes (gradients) in potential, so any constant added or subtracted from the total potential falls out when we take the derivative required to compute the flux.

So how are we going to add elevation to pressure to compute the TMP?  We need to employ common units, and the most natural approach is to use the units of elevation: length. We typically denote the elevation with the variable z which becomes more positive in the upward direction. The pressure potential, however, most naturally has units of force per area(see Table 2.2), and can be calculated within a hydrostatic liquid by multiplying the distance to the water table (say, h) by the liquid density and the acceleration of gravity: \mathrm{\uppsi_{g} = \uprho gh}. Here we have the trick which makes flow calculations work out: we can express the pressure in terms of the equivalent height of a column of water which would generate this pressure. This is the pressure “head,” which we denote h, and has units exactly compatible with elevation. We denote the TMP as H and is equal to the sum of h and z.

To get a feel for how this total mechanical potential works, let’s consider the simplest possible case: a tank of water at rest (Figure 2.X124). Since there is no flow in the tank, the condition which is called “hydrostatic,” we can be sure that the total mechanical potential is the same everywhere (remember water would flow towards lower potential, but we posited that the water was at rest). As mentioned, we may define our datum for pressure and elevation arbitrarily: here we chose the top of the water surface as the zero pressure (denoted with the letter h when in units of head) datum and the bottom of the tank as the zero elevation datum. Clearly elevation increases with height in the tank, while pressure decreases such that the total mechanical, denoted with the letter H when measured in terms of head, potential is constant, consistent with the lack of flow in the tank. Clearly one may convert between units of pressure and equivalent water column height (Table 2.2). A word of caution is in order regarding the assumed density of the fluid which fills the pore space: what about cases of saline solutions, or non-aqueous phase liquids (NAPLs; hydrocarbons etc.)? We must be consistent in our units of pressure measurement, and in this case we will stick to height a column of pure water would be lifted (at 20 deg C) regardless of the fluid under consideration.

Table 2.2 Conversions between common units of pressure. To convert the units given in the left column to that listed on the upper row, multiply by the conversion factor listed.

to:         from: kilo Pascal bar lb./in2 cm water ATM mm Hg
Kilo Pascal 1 0.01 0.145038 10.1972 0.00986923 7.50062
bar 100.000 1 14.5038 1,019.72 0.986923 750.062
lb./in2 6.89476 0.0689747 1 70.3072 0.068046 51.7149
cm water 0.0980665 0.000980665 0.0142233 1 0.000967841 0.735559
ATM 101.325 1.01325 14.6960 1,033.23 1 760
mm Hg 0.133322 0.0013332 0.0193368 1.35951 0.00131579 1

People often wonder about pressures more negative than -1 atmosphere: isn’t this somehow non-physical? We know from our kitchens that boiling occurs when the vapor pressure of the water reaches the atmospheric conditions, so why doesn’t the pore water boil? Don’t feel bad if you came up with these questions; in contemporary literature some authors still wondered if water pressures more negative than -1 ATM were physically possible (Gray and Hassanizadeh, 1991; Miller, 1994; Gray and Hassanizadeh, 1994). However, experiments have shown that liquid water can sustain negative pressures of up to 150 bars before vaporizing at room temperature since the vapor pressure of water is actually much less than one atmosphere at this temperature (e.g., Zheng et al., 1991) supporting Miller’s position. This conclusion is confirmed in trees over 10 m in height, where the pore pressures more negative than -1 bar are essential to the delivery of water to their leaves.  Stable water at such negative pressures only occur in very small bodies of water (e.g., xylem of plants), but at such great negative pressures water will only be found in very small regions, as is discussed below. Thus a few points are clear: (a) negative pressures less than -1 bar arise commonly in porous media due to the forces listed in Table 2.1, and (b) it is not unreasonable to consider the fluid-dynamic behavior of water at pressures more negative than -1 bar.

Figure 2.3a
Figure 2.X124. The terms of the total mechanical potential in a tank of water. The elevation potential in units of head (dotted line) always increases linearly with elevation. Since there is no movement in the tank, we know the total potential (solid line) is the same at all places, and therefore the pressure potential (dashed line) must increase with depth in a manner equal and opposite to the elevation, which agrees with our intuition that pressure goes up linearly with increasing depth in still water.

Example 2-X:  Gravitational Potential

Figure XXX. Self-study of the elevation potential of a soil column to understand the role of selection of datum.  All distance values shown are in mm.

Case 1:

Point A is 100 mm above the datum and point B is 200 mm below the datum. Find the difference in gravitational potential: \mathrm{\upDelta \Psi_{z} = \Psi_{zA} - \Psi_{zB}} between the points.


\mathrm{\Psi_{zA}} = \SI{100}{\mm} and \mathrm{\Psi_{zB}} = \SI{-200}{\mm}, so \mathrm{\Delta \Psi_{z}} = 100 - (-200) = \SI{300}{\mm}.

Case 2:

For the same two points A and B we now change the elevation of the reference to 100 mm above point A.


\mathrm{\Psi_{zA}} = \SI{-100}{\mm} and \mathrm{\Psi_{zB}} = \SI{-400}{\mm}, so \mathrm{\Delta \Psi_{z}} = -100 - (-400) = \SI{300}{\mm}.

Notice we are not here dealing with the possible difference in pressure potential at points A and B, just their elevation potential.  The key point here is that the choice of datum does not change the value of the difference in potential between two points.

Total Soil Water Potential

Definition of Total Soil Water Potential: While the Total Mechanical Potential (TMP, denoted H when measure in units of head) directly gives rise to flow, the energy state of water includes other factors, such as solute concentration and temperature. So the “Total Soil Water Potential” must include all of the quantities that would influence the amount of energy per unit quantity of soil water at equilibrium would be required to take water from some standard reference state to that of the water in the soil. Here we will take the standard reference state to be pure, free water, at the elevation of our datum. The total potential ψT may be expressed as the sum of the contributors to the total energy state:

(16)   \begin{equation*} \psi_{T} = \psi_{H} + \psi_{S} +\psi_{E}+\psi_{Z} + \psi_{P} \end{equation*}

where ψH is the thermal potential, ψs is solute or osmotic potential, ψE is the electrostatic potential, ψz is gravitational potential, and ψp is pressure potential. Outside of engineered systems where electrical potentials are used to drive flow, and at the atomic scale near clay surfaces, the electrostatic potential rarely comes into play, so we will drop it from our accounting for the moment.  For more information on this term we refer you to Revil and Jardani (2013).   The difference in chemical and mechanical potentials between soil water and pure water at the same temperature is known as the soil water potentialw):

(17)   \begin{equation*} \psi_{W} = \psi_{S} + \psi_{P} \end{equation*}

Note that the gravitational and thermal components (ψand ψH) are absent in this definition. Soil water potential is thus the result of inherent properties of soil water itself, and of its physical and chemical interactions with its surroundings, while total potential includes the effects of gravity (an “external” force), or its temperature.  The idea is to define a soil state that can be measured, but where ψW for a sample of soil would not change if you walked from room to room in your laboratory, or put the sample on a high shelf.

Total soil water potential and its components may be expressed in several ways depending on what you choose for your “unit quantity of water.” Potential may be expressed as (i) energy per unit of mass; (ii) energy per unit of volume; or (iii) energy per unit of weight. A summary of the resulting dimensions, common symbols and units are presented in Table 1-1. Note that energy has units of ML2/t2.

Table 1-1: Units, Dimensions and Common Symbols for Potential Energy of Soil Water

Units Symbol Name Dimensions SI Units cgs Units
Energy/Mass μ Chemical Potential L2/t2 J/kg erg/g
Energy/Volume ψ Soil Water Potential, suction, or tension M/(Lt2) N/m2(Pa) erg/cm3
Energy/Weight h Soil Water Head L m cm
(From: Jury et al., Soil Physics, 5th Ed. 1991; reprinted by permission of John Wiley & Sons, Inc.)

Note that only µ has actual units of potential; ψ has units of pressure, and h of head (L) of water. However, the above terminology is widely used in a generic sense in the soil and plant sciences (i.e., potential energy vs. units of potential). The various expressions of soil water energy status are equivalent. Specific relationships between the three potential energy expressions are derived from the relationship between volume and mass of water (\mathrm{M_{w} = \uprho_{w} \cdot V_{w}} where ρw is the density of water), and between mass and weight (\mathrm{W_{w} = g \cdot M_{w}} where g is acceleration of gravity). The relationships may be summarized by:

(18)   \begin{equation*} \mu = \frac{\Psi}{\rho_{w}} = \mathrm{g\ h} \end{equation*}

Example 2.X: Water Potential Conversions

Problem Statement:

Convert a soil water head of -0.01 m to kPa and to J/kg.


(1) Conversion from h=-0.01 m to ψ (kPa):

The conversion equation is \mathrm{\Psi = \uprho_{w}gh} where \sisetup{per-mode = symbol} \mathrm{g} = \SI{9.8}{\meter\per\second\tothe{s}} and \sisetup{per-mode = symbol} \mathrm{\uprho_{w}} = \SI{1000}{\kilogram\per\cubic\meter}.


\mathrm{\uppsi} = -0.01 * 9.8 * \SI{1000}{(\meter * \meter * \kilogram) \per (1\per\square\second * \cubic\meter)}  \rightarrow = \SI{-98}{(\kilogram * \meter\per\square\second) * (1\per\square\meter)}  = \SI{-98}{\pascal}



(note that \sisetup{per-mode = symbol} \si{\pascal = \newton\per\square\meter}, and \sisetup{per-mode = symbol} \si{\newton = kg * \meter\per\square\second}).

(2) Conversion from h=-0.01 m to µ (J/kg):

The conversion equation is \mathrm{\upmu = gh} where \sisetup{per-mode = symbol} \mathrm{g} = \SI{9.8}{\meter\per\square\second}, hence

\sisetup{per-mode = symbol} \mathrm{\upmu} = -0.01 * \SI{9.8}{(\meter * \meter\per\square\second)} \rightarrow = \SI{-0.098}{(\meter * \meter\per\square\second});

we multiply and divide the expression by

\sisetup{per-mode = symbol} \si{kg} \rightarrow = \SI{-0.098}{\kilogram * \meter * \meter\per(\kilogram * \square\second)}

which is

\sisetup{per-mode = symbol} \SI{-0.098}{\joule\per\kilogram}

(note that a Joule = \sisetup{per-mode = symbol} \si{\newton * \meter = (\kilogram * \meter\per\square\second) * \meter}).

When do we employ total potential, versus mechanical potential, in the estimation of fluxes from high to low potential? As we discussed above, the osmotic potential can cause movement of water when there is an interface in the system which constrains the flux of selected constituents. For example, if there were two side-by-side cups of water solutions in a sealed box which could exchange water only in the vapor phase, then the air-water interfaces would act as constraints on the non-volatile solutes, while letting water be freely exchanged.  If left long enough the water would move in the vapor phase until the water came to have the same total potential in the two cups. For instance, if one was pure water, and the other had a dissolved salt, then the pure water would tend to lose water in the vapor phase which would then condense into the saline solution, diluting the salt, and bringing the total potential of the two vessels closer to equal. Another important example is that of plant turgor: in contact with pure water, plants can obtain high cell pressures by having high dissolved concentrations in their tissue (low potential), so that the equilibrium total potential is obtained when the internal pressure in the tissue is high. Similarly, if the salt concentration in soil is too high, plants cannot create an energetic gradient to drive the saline water into the plant, and thus can wilt.

You will often hear of water described to be under “tension,” “matric potential” or “suction.” These terms all refer to the same quantity, the strength of attraction of water to a medium.This term is not rigorously tied to a particular physical mechanism, though some authors employ this term to identify the energy of interaction between water and solid surfaces (recall Figure 2.X illustrating the high energy of this interaction by way of the temperatures exceeding 700 C to remove all of this water – for a classic treatment of this see Brown and Montgomery, 1913, “Dehydration of Clays,” Bureau of Standards Technologic Paper no. 21, Washington, D.C.: https://nvlpubs.nist.gov/nistpubs/nbstechnologic/nbstechnologicpaperT21.pdf). While this is a perfectly fair energy term to quantify, it is important to recognize that, in general, solid surfaces will be at some intermediate state of water absorption, so any theoretical maximum value is rarely relevant. Generally speaking the energy per unit area of this interface is far in excess of the air-water surface tension, which theoretically should lead water to spread across all such surfaces.

Hydrostatics in Porous Media

Having defined the parameterization of the porous media, and the units of the liquid phase potential, we can now discuss the interactions between soil and fluids. In principle there could be any number of fluids in these spaces, say m, each of them filling some partial volume of the pore space; V1, V2, …Vm. The only constraint on this system is that

(19)   \begin{equation*} \sum\limits_{i = 1}^m \mathrm{V_{i} = n} \end{equation*}

where we recall that n is the porosity. Fluids of primary interest are air and water; however other gases and non-aqueous phase liquids (e.g., gasoline and diesel in both liquid and gaseous forms) are also of significant interest. The focus of this text are on water and air, but the formulations are sufficiently general to allow direct extension to other fluids.

Why doesn’t all the water just fall out of the porous media? There are at least four groups of forces contributing to the retention of water in soil, which are enumerated here in order of decreasing strength. Each of these forces has a connection to the structure of water, which has a polar nature, slightly negatively charged on the side where oxygen projects, and positively charged on the side where hydrogen is prominent.  This means that water molecules have a highly dipolar electrostatic structure.  Water molecules then develop a hydrogen bond between neighbors between their oxygen and hydrogen.  Among other things, this leads to high bonding energy between water molecules, which leads to high surface tension, high latent heat of vaporization, and high solubility of ions.  To summarize:

  1. The hydrogen of water is attracted to the negative surface charge of mineral surfaces (Van der Waals attraction).
  2. The periodic structure of the clay surfaces gives rise to an electrostatic dipole which results in an attractive force to the water dipole.
  3. Osmotic force, caused by ionic concentration near charged surfaces, hold water.
  4. Surface tension at water/air/solid interfaces maintains macroscopic units of water in pore spaces.

Table 2.1 Length scale (r) of forces contributing to retention of water in soils.

Force Attraction law Range of influence
Van der Waals 1/r7 1 Å
Electrostatic – Dipole 1/r2 – ion

1/r6 – dipole

5 Å
Osmotic (double layer) 1/r 100 Å
Surface Tension 1/r 1 mm

These forces are discussed in detail in Sposito (1989) and Kirkham and Powers (1972) for readers who would like to learn more about these fascinating and fundamental relations. Attention in this text is focused on transport phenomena in the vadose zone. The first three forces listed in Table 2.1 act only at very short range, and essentially immobilize water. In contrast, surface tension produces effects which affect water in bulk, and is therefore more influential in transport.

Surface Tension

Surface tension is a measure of the energy required to generate an interface between two phases.  Here we will be most concerned with the air-water interface, as shown in Figure 2.X2, but the concept applies to solid-gas interfaces, and liquid-solid interfaces.  An intuitive understanding of surface tension can be obtained from a simple thought experiment. Imagine a block of water in a container that can be split in two (Figure 2.X3). For convenience think of each water molecule as being organized in a cubic lattice structure, with six bonds. Now we quickly split this block of water into two halves. The molecules on the new air/water surfaces are bound to only five of their neighbors. It took energy to break these bonds, which will be manifest as a free surface energy. Since the water surface has a constant number of molecules on its surface per unit area, the energy required to create these surfaces is directly related to the surface area created. This would suggest that surface tension should have units of energy per unit area, or force per unit length.

Figure 2.X2.  Illustration of the tightly organized water molecules in the liquid state, and the constant adsorption of water vapor and desorption of liquid water from the air-water interface. The actual structure of water at the interface differs considerably from either the pure liquid or gaseous phases (e.g., Tong, Y., T. Kampfratha and R.K. Campen, Experimentally probing the libration of interfacial water: the rotational potential of water is stiffer at the air/water interface than in bulk liquid, Phys. Chem. Chem. Phys., 2016,18, 18424-18430, DOI 10.1039/C6CP01004K) in that there is a transition zone where the density of water moves between the liquid and vapor values.


Figure 2.3
Figure 2.X3. Conceptual model of potential energy of broken bonds at liquid surface. The surface area of system (a) is 10 square units, while that of system (b) is 12 square units. System (b) has greater potential energy due to the broken surface bonds, as illustrated for an open vessel of water in figure (c).


It is easy to measure surface tension using the sliding wire set-up shown in Figure 2.X4. If the force measured on the scale is F, and the width of the film is L, then the surface tension, σ, is simply

(20)   \begin{equation*} \mathrm{\upsigma = \frac{F} {2L}} \end{equation*}

Figure 2.4
Figure 2.X4: Wire Frame set-up for measuring surface tension of air-liquid interfaces.

Wait a minute! How did that factor of 2 sneak into Equation 20? Simple: since the film is open on two sides, we are creating two air/water interfaces, and thus the effective length of a single-sided free surface is 2L. Note that you could also keep the entire problem in terms of energy per unit area in keeping with the thought experiment of the previous paragraph: carefully record the product of force required to stretch the film multiplied by the distance over which that force is exerted (to obtain the energy input) and divide by the area (l x L). Table 2.3 lists some typical values of surface tension for fluid pairs.

Table 2.3 Some values of surface tensions for various interfaces (taken from the CRC Handbook of Chem. and Physics 68th ed., 1987).

Substance Formula Contacting Temp (℃) σ


Water H2O air 20  72.75
Acetone C3H6O air 20 23.7
Benzene C6H6  







Carbon Tetrachloride CCl4 vapor






Mercury Hg water 20 375

Classical Derivation of Laplace’s Equation

We are concerned with how surface tension manifests itself in a porous medium. In particular, we would like to know the static fluid pressures which result from surface tension acting on fluid surfaces, which are curved within a porous medium. Consider the infinitesimal curved fluid surface which is shown in Figure 2.5. For the sake of generality, we allow for curvature along two axes of radius r1 and r2. To determine the pressure p acting on the surface we simply add up the forces acting on the interface in the vertical direction. Since the system is at equilibrium, these forces must add to zero. The forces acting upward (which we will take as positive) result from p2 and r2

(21)   \begin{equation*} \mathrm{F_{up} = p_{2}\updelta S_{1} \updelta S_{2} + 2 \upsigma \updelta S_{1} \sin \left(\frac{\updelta \uptheta_{2}} {2}\right)} \end{equation*}

where σ is the surface tension of the interface. Now since this is a very small patch, the angle δθ2 is very small, and we may use the small angle approximation for sine:

(22)   \begin{equation*} \sin \left(\frac{\updelta \uptheta_{2}} {2}\right) \approx \frac{\updelta \uptheta_{2}} {2} \approx \frac{\updelta \text{S}_{2}} {2\text{r}_{2}} \end{equation*}

By symmetry, the expression for Fdown is identical to Equation 21, substituting appropriate subscripts and signs. The net force on the surface is then given by

(23)   \begin{align*} \mathrm{F_{net}} = &0 = \text{upward forces} - \text{downward forces} \nonumber \\ \nonumber \\ &= \mathrm{p_{2} \updelta S_{1} \updelta S_{2} + 2 \upsigma \updelta S_{1}\frac{\updelta S_{2}} {2r_{2}} - p_{1} \updelta S_{2} \updelta S_{1} - 2 \upsigma \updelta S_{2} \frac{\updelta S_{1}} {2r_{1}}} \end{align*}

which, canceling δθ1δθ2 from all terms, simplifies to

(24)   \begin{equation*} \mathrm{p_{2} - p_{1} - \frac{\upsigma} {r_{1}} + \frac{\upsigma} {r_{2}} = 0} \end{equation*}


(25)   \begin{equation*} \mathrm{p_{2} - p_{1} = \upsigma \left(\frac{1} {r_{1}} - \frac{1} {r_{2}}\right)} \end{equation*}

which is Laplace’s Equation for pressure across an interface at a saddle point (an anticlastic surface). Notice that if both radii were of the same sign and magnitude (i.e., the surface was spherical) we have r1 = – r2 = R in Figure 2.5, and Equation 25 would reduce to

(26)   \begin{equation*} \mathrm{\upDelta p = \frac{2\upsigma} {R}} \end{equation*}

which is the familiar form of Laplace’s equation for pressure across a spherical fluid interface. Of course this is the same expression one obtains for the pressure behind the gas/liquid interface in a capillary tube of radius R. For the general case where \mathrm{r_{1} \neq r_{2}}, we define the reduced radius of curvature, R, as

Figure 2.5
Figure 2.5. Definition sketches for a doubly curved (anticlastic) liquid surface.

(27)   \begin{equation*} \mathrm{\frac{2} {R} = \left(\frac{1} {r_{1}} + \frac{1} {r_{2}}\right)} \end{equation*}


(28)   \begin{equation*} \mathrm{R = \frac{2r_{1} r_{2}} {r_{1} + r_{2}}} \end{equation*}

where we have used a sign convention for r1 and r2 such that they are taken to be positive if measured in the wetting fluid and negative if measured in the non-wetting fluid (typically air).  Note that an element of water in a porous media often displays complex curved interfaces.  Even the simplest case of two spheres in contact with water between has two characteristic radii, as shown in Figure 2.X

(a) Illustration of the combined convex and concave gas-liquid interface at the contact of two spherical particles as well as the progression of the water shape as a function of pressure (b) (from XXX with permission).

Introduction to the of Fluid Interfaces

The derivation of Laplace’s equation provided in section is completely acceptable, but in this case there is good reason to look at this problem from another perspective, that of thermodynamics, to gain further insight into the workings of fluid interfaces. In the following sections we shall re-derive Laplace’s equation, using a somewhat less intuitive and more formal approach. This will be a warm-up, for then we will push thermodynamics a bit further to obtain some very handy additional results.

The Thermodynamics of Flat Interfaces

Consider a small, partially filled container of water with a flat water-air interface. The system has water in three phases: liquid, gaseous and as a taut interfacial surface, denoted by subscripts ‘l’, ‘g’, and ‘a’, respectively. We will denote the surface area as Σ and the surface tension at the air/water interface as σ. Since there is no barrier to water moving between phases, the chemical potential of the water in each phase is the same, i.e.,

(29)   \begin{equation*} \mathrm{\upmu_{g} = \upmu_{a} = \upmu_{l} = \upmu_{\infty}} \end{equation*}

The subscript ∞ is used to indicate values of intensive parameters which are constant in all three phases (recall intensive parameters are those which do not depend upon the amount of material in question, in contrast to extensive parameters which do depend on the amount, or ‘extent’ of the material). We will consider a system at constant temperature, T. We will use Ni to denote the number of atoms, Si the entropy, Ei the internal energy, Pi the pressure and Vi the volume in phase i.

Considering the gas and liquid phases, we have three fundamental differential forms

(30)   \begin{equation*} \mathrm{TdS_{g} = dE_{g} + P_{\infty} dV_{g} - \upmu_{\infty} dN_{g}} \qquad \qquad \qquad \text{(gas)} \end{equation*}

(31)   \begin{equation*} \mathrm{TdS_{1} = dE_{1} + P_{\infty} dV_{1} - \upmu_{\infty} dN_{1} - \upsigma d\Sigma} \qquad \qquad \qquad \text{(liquid)} \end{equation*}

(32)   \begin{equation*} \mathrm{TdS_{a} = dE_{a} - \upsigma d \Sigma} \qquad \qquad \qquad \text{(surface)} \end{equation*}

The total energy and entropy for the system are simply the sum of the components, i.e.

(33)   \begin{equation*} \mathrm{S = S_{a} + S_{g} + S_{1}} \end{equation*}

(34)   \begin{equation*} \mathrm{E = E_{a} + E_{g} + E_{1}} \end{equation*}

Since the surface is two-dimensional, we take the number of atoms in the surface to be zero in comparison to the atoms in the three-dimensional volumes of gas and liquid. So we may write

(35)   \begin{equation*} \mathrm{N=N_{1}+N_{g}} \end{equation*}

Thus the fundamental differential form for the system is given by

(36)   \begin{equation*} \mathrm{TdS = dE + P_{\infty}dV - \upmu_{\infty} dN - \upsigma d \Sigma} \qquad \text{(system)} \end{equation*}

Employing the Gibbs-Duhem relationship (i.e., for an exact differential, the differentiation may be shifted from the extensive to intensive variables maintaining equality) we find

(37)   \begin{equation*} \mathrm{S_{a} dT = \Sigma d \upsigma} \end{equation*}


(38)   \begin{equation*} \mathrm{\frac{S_{a}} {\Sigma} = \frac{d\upsigma} {dT}} \end{equation*}

which is the equation of state for the surface phase (analogous to the equation of state for an ideal gas \mathrm{PV = nRT}), and relates the temperature dependence of surface tension to the magnitude of the entropy of the surface.

Laplace’s Equation from a sphere in Space

Now let’s consider the effect of a curved air-water interface. Now Pg and Pl are not equal, while μg and μl remain equal. Considering a droplet floating in gas filled space of constant volume, the fundamental differential form for the system is given by

(39)   \begin{equation*} \mathrm{TdS = dE + P_{g} dVg + P_{1} dV_{1} - \upmu_{g} dN_{g} - \upmu_{1} dN_{1} - \upsigma d \Sigma} \end{equation*}

Considering a very small transfer of volume, δV, between the two phases, we see the chemical potential terms would be equal and opposite, that the total change in energy in the system is unchanged (we are doing no work on the system), and the entropy would be constant to first order. Since the volume of the system is held constant from 39 we find

(40)   \begin{equation*} \mathrm{\left(P_{1} - P_{g}\right) \updelta V - \upsigma \updelta \Sigma = 0} \end{equation*}

or defining Pd as the pressure difference across the interface of the droplet, we have

(41)   \begin{equation*} \mathrm{P_{d} = \upsigma} \frac{\delta \mathrm{S}} {\delta \mathrm{V}} \end{equation*}

We can calculate the differential on the right-hand side by noting that for a sphere \mathrm{V = (4\uppi r^{3}/3)} and \mathrm{\Sigma = 4\uppi r^{2}} , so \mathrm{\updelta \Sigma / \updelta V = 2/r}, hence

(42)   \begin{equation*} \mathrm{P_{d} = \frac{2\upsigma} {r}} \end{equation*}

which is Laplace’s equation for the pressure across a curved interface where the two characteristic radii are equal (see 26).

Although you probably have had it with this equation by now, it is only fair to show the simple way to obtain this result. Consider the pressure balance across the middle of the droplet. We know the surface tension of the water about the center of the droplet must equal the pressure exerted across the area of the droplet by the liquid. The area of the droplet at its midpoint is πr2 at pressure Pd, while the length of surface applying this pressure is 2πr at tension σ. Equating these we obtain

(43)   \begin{align*} &\mathrm{P_{d} \uppi r^{2} = 2 \uppi r \upsigma} \\ \nonumber \\ &\mathrm{P_{d} = \frac{2 \upsigma} {r}} \end{align*}

as expected.

Across Curved Interfaces

The curved interface also affects the vapor pressure of the liquid species. Suppose we have our spherical water droplet in a fixed volume at equilibrium with the vapor phase. The chemical potential of the water in these two phases will be equal, and will remain equal through any reversible process, i.e.,

(44)   \begin{equation*} \mathrm{\upmu_{1} = \upmu_{g}} \end{equation*}


(45)   \begin{equation*} \mathrm{\updelta \upmu_{1} = \updelta \upmu_{g}} \end{equation*}

Writing the fundamental differential forms for the gas and liquid phases for a system with one mole of each species (following equations 30 and 31) we have

(46)   \begin{equation*} \mathrm{TdS_{g} = dE_{g} + P_{g} dV_{g} - \upmu_{g} dN_{g}} \qquad \qquad \qquad \text{(gas)} \end{equation*}

(47)   \begin{equation*} \mathrm{TdS_{1} = dE_{1} + P_{1} dV_{1} - \upmu_{1} dN_{1}} \qquad \qquad \qquad \text{(liquid)} \end{equation*}

As before, the Gibbs-Duhem relation of the fundamental differential forms follows the principle that for an exact differential form, equality is maintained if differentiation for each term is exchanged between intensive and extensive variables, hence 46 and ?? may be written

(48)   \begin{equation*} \mathrm{S_{g} dT = V_{g} dP_{g} - N_{g} d \upmu_{g}} \qquad \qquad \qquad \text{(gas)} \end{equation*}

(49)   \begin{equation*} \mathrm{S_{1} dT = V_{1} dP_{1} - N_{1} d \upmu_{1}} \qquad \qquad \qquad \text{(liquid)} \end{equation*}

Dividing 48 by Ng and 49 by Nl we get the equations in a per mole basis. By further assuming that temperature and surface tension remain constant, we obtain

(50)   \begin{equation*} \mathrm{v_{g} dP_{g} = d \upmu_{g}} \qquad \qquad \qquad \text{(gas)} \end{equation*}

(51)   \begin{equation*} \mathrm{v_{1} dP_{1} = d \upmu_{1}} \qquad \qquad \qquad \text{(liquid)} \end{equation*}

where v in lower case indicates the volume per mole. Now employing ?? we find

(52)   \begin{equation*} \mathrm{v_{g} dP_{g} = v_{1} dP_{1}} \end{equation*}

which may be written

(53)   \begin{equation*} \mathrm{d \left(P_{1} - P_{g} \right) = \frac{v_{g} - v_{1}} {v_{1}} dP_{g}} \end{equation*}

Using Laplace’s equation (e.g., 43) this may be re-written as

(54)   \begin{equation*} \mathrm{d \left(\frac{2 \upsigma} {r} \right) = \frac{v_{g} - v_{1}} {v_{1}} dP_{g}} \end{equation*}


(55)   \begin{equation*} \mathrm{dP_{d} = \frac{v_{g} - v_{1}} {v_{1}} dP_{g}} \end{equation*}

We note that vl is approximately four orders of magnitude less than vg, so we will suppose that \mathrm{(v_{g} - v_{l}) / v_{l}} can be accurately approximated by \mathrm{v_{g} / v_{l}}, and we obtain

(56)   \begin{equation*} \mathrm{d \left(\frac{2 \upsigma} {r}\right) = \frac{v_{g}} {v_{1}} dP_{g}} \end{equation*}

Now assuming the vapor phase acts like an ideal gas, \mathrm{P_{g}v_{g} = RT}, whence 56 becomes

(57)   \begin{equation*} \mathrm{d \left(\frac{2 \upsigma} {r}\right) = \frac{RT} {P_{g} v_{1}} dP_{g}} \end{equation*}

This may be integrated from a flat interface (\mathrm{ r = \infty}) to that with radius \mathrm{r} to obtain

(58)   \begin{equation*} \mathrm{\frac{2 \upsigma} {r} = \frac{RT} {v_{1}} \ln \left(\frac{P_{g}} {P_{\infty}}\right)} \end{equation*}

where P is the vapor pressure of water at temperature T (Pg/P is the relative humidity of the soil air system). Taking R̄ as the specific gas constant for water (i.e., R̄ = R/vl), and noting that the left-hand side of 58 is just Pd, this may also be written

(59)   \begin{equation*} \mathrm{P_{d} = \bar{R} T \ln \left(\frac{P_{g}} {P_{\infty}}\right)} \end{equation*}

where Pd is the pressure potential of the water. This very useful relation is known as the psychrometric equation, which provides the relationship between the humidity of the air in equilibrium with an air-water interface and the potential of the water.  This allows the determination of very high negative potentials through measurement of the vapor pressure (or relative humidity) of water in porous media. For instance, at a pressure of -15 Bar (the permanent wilting point of many plants), Pg/P is 0.99. Measurement of Pg/P may be made in a number of ways. Most commonly a small thermocouple is cooled by application of a reverse current while its temperature is read with a second thermocouple. As the temperature of the first thermocouple attains the dew point of the vapor, its temperature decline is sharply reduced due to the energy of condensation of water on the device. Knowing the dew point, it is straightforward to obtain the relative humidity with good precision (Rawlins and Campbell, 1986).

Temperature Dependence of Surface Tension

It is often overlooked that all the measurements we take regarding water/media interactions are strongly temperature dependent. In many of the equations derived above we have treated surface tension as a constant property of a fluid. This is certainly not the case. The temperature dependence of air-water surface tension is illustrated below in Figure 2.6 where we see that for this pair, surface tension decreases at approximately one percent per 4℃, which could give rise to discrepancies of greater than 5% between common laboratory and field conditions, or even through diurnal fluctuations.

Figure 2.6
Figure 2.6. The temperature dependence of air-water surface tension (Data from the 68th ed. of the Handbook of Chemistry and Physics, 1987).

The Characteristic Curves

Now for the real excitement: putting the solid and fluid together! To do this, we will undertake yet another simple thought experiment. Take a funnel, and seal a flat porous plate in it to make a Buchner funnel, as shown in Figure 2.7. Select the porous plate so that its pores, of radius r, are much smaller than the pores in the soil sample. Now attach a tube to the bottom of the funnel, and fill the region below the porous plate to the end of the tube with water. This tube is fed from a beaker on a precision balance, which keeps track of the amount of water in the soil sample. By wetting the porous plate, we know from Equation 26 that air cannot enter the porous plate until a suction greater than Pentry

(60)   \begin{equation*} \mathrm{P_{entry} = \frac{2 \upsigma} {r_{max}}} \end{equation*}

is applied to the plate. Since the porous plate has pore radii smaller than the pore in the soil sample, we know that the plate will remain saturated as we draw the water from the sample. We now place a thin piece of dry soil on the porous funnel surface. The experiment will take place over three stages: (1) wetting, (2) drying, and (3) re-wetting. For stage (1), we incrementally bring the end of the tube up, finally attaining the elevation of the soil. At each lift we measure the amount of water drawn up by the soil, with pressures ranging from Pentry to 0. Each measurement is taken allowing the system to come to a steady state. Measuring elevation, H, as positive upward, the pressure applied to the water in this soil will be given by

(61)   \begin{equation*} \mathrm{P_{soil} = \uprho_{w} gh}} \end{equation*}

where ρw is the density of water and h is the pressure head of the water in the sample, which we will assume is much larger than the thickness of the sample l. In stage (2) the end of the tube is incrementally lowered to apply a negative head to the water in the soil. Finally, in stage (3) we again lift the beaker to an elevation slightly below the height of the soil sample keeping track of the amount of water taken up by the sample. A typical plot illustrating how the soil sample’s volumetric water content (equation \ref{eq:volumetric_liquid}) changes in response to pressure changes for this type of experiment is shown in Figure 2.8, with a few more lifts and drops thrown in for good measure. The vertical axis in Figure 2.8 represents the pressure in the sample, as reflected by the hanging height of water, h. The horizontal axis shows the soil’s volumetric moisture content.

The equilibrium relationship between pressure and moisture content is called the soil water characteristic (SWC). The SWC is an important soil property related to the distribution of pore space (sizes, interconnectedness), which is strongly affected by texture and structure, as well as related factors including organic matter content. The SWC is a primary hydraulic property required for modeling water flow, for irrigation management, and for many additional applications related to managing or predicting water behavior in the porous system. A SWC is a highly nonlinear function and is relatively difficult to obtain accurately. Because the pressure potential extends over several orders of magnitudes for the range of water contents encountered in the field, the potential is often plotted on a logarithmic scale (Figure 2.X).

Figure 2.X SWC curves for soils of different textures demonstrating the reduced porosity (saturated water content) and the relatively large pores for coarse-textured media such as sand, as well as the extreme range of potentials required to drain soils of different texture.

This simple experiment illustrates most of the physics which control the retention and movement of fluids through porous media, so it is worth considering these results carefully. First consider the wetting curve alone, spanning the points labeled (1), (2) and (3) in Figure (2.8). We see that more water is taken up by the soil (i.e., the moisture content increases) as the beaker comes closer to the elevation of the soil (i.e., as the pressure of the feed water increases).

One can gain a basic understanding of the characteristic S-shape of this curve by relating it to the soil particle and pore size distributions. Consider a soil made up of particles with size distributed between dmin and dmax as shown in Figure 2.9. By passing the sample through successive sieves, from very fine to coarse, we may evaluate the percentage of the grains which have size greater than any given size, as plotted in Figure (2.9a) (see section 1.3.2). Without stretching your imagination too far, it is plausible that the pore size distribution in the original soil is quite similar to the particle size distribution. Now the ordinate goes from mass of particles, to volume of pores. The plot of pore size distribution which is obtained is quite similar to that of particle size, as shown in Figure (2.9b).

This is often conceptualized in terms of the soil being made of many differing pore sizes classes which could be represented as equivalent parallel capillary tubes of many radii (Figure 2.X).

Figure 2.X Capillary tube conceptual models for the pores in soils. (a) illustrates the notion that the pores of larger diameter draw water to lesser heights, while fine pores draw water higher. Thus, a coarse soil will draw in water with less pressure than will a fine soil. (b) considers a single vertical soil sample at equilibrium with water provided to the lower boundary, with a variety of pore sizes represented by capillary tubes. This shows that all pores will be filled at 0 pressure, then as the the pressure potential becomes more negative with increasing elevation (recall that at equilibrium the total potential is constant in the sample), only ever smaller pores can retain water.
(a) (b)

Using Equation 26, we can equate pore size to the capillary filling pressure of each pore, and re-plot the pore size distribution as a filling pressure plot; Figure (2.9c). Notice the volume of pores which are filled in Figure (2.9c) is simply the degree of saturation. At last, in Figure (2.9d), we have obtained the wetting curve from the particle size distribution. The idea here is to point out the connection between the wetting characteristic curve and the physical properties of the medium. A similar argument is often used when trying to determine how a medium will drain using the particle size distribution (e.g., Arya and Paris, 1981).

Figure 2.7
Figure 2.7. (a) Conceptual illustration of Bukner funnel experiment to determine characteristic draining curve; (b) close-up view of idealized porous plate.
Figure 2.8
Figure 2.8. The main characteristic curves for wetting, draining, and one set of scanning curves.
Figure 2.9
Figure 2.9. Conceptual relationship between (a) particle size distribution; (b) pore size distribution; (c) pore volume occupied vs pore filling pressure; and (d) the characteristic wetting curve.

Now let’s follow the process of draining, going from points 3 to 4 and 5 in Figure (2.8). Why isn’t this curve following right along the wetting curve? To understand this we shall consider the filling and emptying of a single pore, as shown in Figure (2.10). The idea represented in this figure again follows directly from Equation 26, and were introduced by Haines in 1930.

The wetting process for this pore is illustrated in stages 1 through 3 of Figure 2.8 and (b) of Figure 2.10, while the draining process is illustrated in stages 4 through 7 of Figure 2.8 and (c) and (d) of Figure 2.10. For any given pressure, there are two states that this pore may be in, depending on whether it is in the process of filling or draining. In this sense the saturation of the pore is a function of both the pressure and the history of the pore. For this reason, the wetting and drying curves for soil are referred to as “hysteretic.

Figure 2.10 shows that the filling process is held up until the pressure at the soil is sufficient to fill the largest pores in the medium. On the other hand, once the medium is filled, it cannot drain significantly until the pressure becomes sufficiently negative to pull the water out of the pore necks that lead into the medium. Clearly between the pressures of -2σ/r2 and -2σ/r1 the moisture content is a function of not only the pore pressure, but what the history of wetting has been.

Based on the ideas contained in Figure (2.9), hw can be identified as the pressure at which the largest group of pore bodies fill. From Figure (2.9) and Equation 26, this indicates that the largest pores in the medium are of radius rmax = 2σ/hw. Similarly, the pressure ha is dictated by the size of the typical pore throats, which have a radius of approximately 2σ/ha.

Why doesn’t the soil drain completely, and have point #5 in Figure (2.8) return to the dry state? There are two phenomena which account for this behavior, with relative magnitude depending on the chemical and physical properties of the medium. First, any water which entered into chemical bonds with charged surfaces in the medium cannot be removed through simple fluid tension: the bonds holding them are as strong as surface tension, and so they cannot be ripped away by the bulk fluid. Furthermore, even in an inert porous medium, fluid can be held in the very small radius regions at particle contacts. These isolated donut shaped pockets of water, referred to as “pendular,” remain in place, no longer hydraulically connected to the fluid reservoir (as shown in Figure 2.10 (d)).

Figure 2.10
Figure 2.10. Illustration of a “Haines Jump” cycle for a set of equal size pores with equal size necks connecting them. In (1) we see that no water enters until the head becomes greater than -2σ/r1. When this pressure is exceeded, the pores will suddenly fill (2). In the draining process ((4), (5), and (6)) first the outer pores will drain, and finally, when the head becomes less than -2σ/r2 all but the isolated pores (shown as the horizontal ones here) will suddenly drain as air can finally enter the necks.

Contact Angle; Sources and Effects

The main hysteretic mechanism illustrated in Figure 2.10 is referred to as the “ink bottle effect,” due to the similarity of the behavior of bottles of ink with very narrow spouts. Another source of hysteretic behavior is hysteresis in the contact angle. Although contact angle was mentioned before in passing, it is now time to take a closer look at this interesting corner of fluid behavior.

Figure 2.11 illustrates a cross-section of a straight-line contact between a fluid and a solid. What should the angle of contact between the solid and fluid be, and why? Well, since the system is at equilibrium, we should be able to balance the forces at the point of contact. Considering only the component of forces acting in the horizontal direction, we find

Figure 2.11
Figure 2.11. (a) A Cross-section of a line contact between a solid plane with a line of fluid/air contact at equilibrium. (b) The force balance at the point of contact.

(62)   \begin{equation*} \mathrm{\Sigma F = 0 = F_{sl} - F_{sg} + F_{lg} \cos \upgamma} \end{equation*}

where Fsg is the solid-gas surface force per unit length, Fsl is the solid-liquid force, Flg is the liquid-gas surface force and γ is the contact angle. Per unit length, Fij = σij, so this may be put in terms of the relative surface tensions

(63)   \begin{equation*} \mathrm{\upsigma_{sg} = \upsigma_{sl} + \upsigma_{lg} \cos \upgamma} \end{equation*}

which is the Youngs-Laplace equation; the relationship between the contact angle and the three surface tensions. Solving [2.56] for the contact angle we find

(64)   \begin{equation*} \mathrm{\upgamma = cos^{-1} \left(\frac{\upsigma_{sg} - \upsigma_{sl}} {\upsigma_{lg}} \right)} \end{equation*}

Care needs to be taken when applying these results, as there are some physical limits on possible values of γ: the contact angle is bounded by 0\degree \leq \gamma \leq 180\degree. So if the operand of cos-1 is greater than 1, then g will be 0 o, while if the value is less than -1, the value will be 180o. Since it is often true that \mathrm{(\upsigma_{sl} - \upsigma_{sg}) > \upsigma_{lg}} for water, the contact angle for water going into geologic material is often taken to be 0 o.

Figure 2.XX: Illustration of wetting and non-wetting surfaces. (a) Wettable silt soil surface (γ ~ 0o) (b) Treated water-repellant silt soil (γ ~ 70o) (Bachmann et al., 2000) ]

Now how does this contact angle relate to hysteresis? As shown in Figure 2.12, the contact angle often exhibits a different value in the advancing and receding cases. This is the reason a drop of water on a flat plate will not start to move as soon as you tilt the plate: more energy is required to remove the water from the trailing edge of the plate than is given up by the sum of the gravitational potential plus the energy released wetting the plate. This is sometimes referred to as the rain-drop effect, perhaps referring to how a raindrop moves down a window pane. The relevance of the rain-drop effect to capillary hysteresis is simply an extension of the observation regarding the plate and drop: a medium will retain water more vigorously than it will absorb water (Bartell and Wooley, 1933).

That this rain-drop effect exists is well known to any child who has bemoaned a rainy day at the window (without it all the droplets would run off the pane), but what is the physical basis? Consider the effect which might result from a scratch on a solid surface as shown in Figure 2.13. At the microscopic level the Youngs-Laplace equation is adhered to, but from a macroscopic point of view, the drip cannot advance until the apparent contact angle is quite large. Similarly, upon retreat the macroscopic contact angle will be much smaller than the true microscopic magnitude (Figure 2.13). The same pair of effects would result from a smudge of oil or other contamination which caused a local reduction in the microscopic contact angle.



Figure 2.12
Figure 2.12. The rain-drop effect of contact angle hysteresis between advancing and receding contact angle, geometrically (left) and as seen (right).








Figure 2.13
Figure 2.13. Illustration of the effect of a microscopic defect on the apparent contact angle for advancing and receding fluid. Although the microscopic contact angle is constant (γ1), the advancing macroscopic contact angle is γ2 and the receding wetting front with macroscopic contact angle γ3 as a result of the surface defect.

A fine cylindrical tube (capillary tube) is often taken as the simplest model for capillary attraction. Given its prominence in describing capillary processes, and because it is a useful conceptual model, we now discuss this system. Figure 2.14 shows two possible capillary tubes. The tube on the left has a zero-contact angle between the tube and the liquid, while the tube on the right has a non-zero contact angle.

Figure 2.14. Examples of capillary tubes or radius r drawing up liquid from a reservoir. (a) With the same liquid, but capillaries with different contact angle (the tube on the left has a zero contact angle between the liquid and tube and the right tube is non-zero) we see the capillary rise is reduced; (b) The maximum capillary rise is dictated by the relative surface energy of the air-solid and liquid-solid interfaces. If the gas-liquid interface had infinite energy, the interface in the tube would be a straight line with a 90-degree contact angle and the interface at the height dictated by the solid interfacial energies. If it happened that the solid-gas energy was very high, the liquid would coat the entire inside of the tube, and the visible gas-liquid interface would have a zero contact angle to the solid. Intermediate values of the gas-liquid interfacial energy give rise to contact angles between 0 and 90 degrees. The lowest interface shown in (b) illustrates a reduction in solid-gas interfacial energy (which reduced the capillary rise), but yet with a high gas-liquid interfacial energy, meaning the contact angle is between 0 and 90 degrees.

So why does the water go up into the tubes? Because energy is liberated when this liquid coats solid surfaces. The ultimate height of the capillary surface is that which balances the energy required to pull the liquid up against gravity and the energy liberated when the liquid coats the solid surfaces. Thus, the derivative of the gravitational energy and surface tension terms are equal. The gravitational (EG) and surface (ES) energies of the system are

(65)   \begin{equation*} \mathrm{E_{G} = \int^{h}_{0} \uprho g \uppi r^{2} h'dh' = \frac{\uprho g \uppi r^{2} h^{2}} {2}} \end{equation*}


(66)   \begin{equation*} \mathrm{E_{S} = \int^{h}_{0} 2 \uppi r \left(\upsigma_{sl} - \upsigma_{sg} \right) dh' = 2 \uppi r \left(\upsigma_{sl} - \upsigma_{sg} \right) h} \end{equation*}

taking the reservoir surface as a datum, for liquid of density ρ. A stable elevation demands that the rate of change in these equal and opposite with change in height, so we find

    \begin{equation*} \mathrm{\uprho g \uppi r^{2} h = -2 \uppi r \left(\upsigma_{sl} - \upsigma_{sg} \right)} \end{equation*}


(67)   \begin{equation*} \mathrm{h = -\frac{2 \left(\upsigma_{sl} - \upsigma_{sg}\right)} {\uprho gr}} \end{equation*}

Which is quite akin to our good friend the Laplace equation, but only the same if . The tricky part here, which is often unappreciated, is that the gas-liquid interfacial energy does not show up in this result, nor does the contact angle. Many people write 67 in terms of the gas-liquid interfacial energy and the contact angle using ?? but this really obscures the important physical nature of the process, as pointed out by Selker and Schroth (1998).

We should look more closely at the case of the zero contact angle, which contains two special cases: (1) \mathrm{\upsigma - \upsigma_{sg} = \upsigma_{lg}} and (2) \mathrm{\upsigma - \upsigma_{sg} = \upsigma_{lg}}. In case (1), the result of 67 is applicable. But in the case (2), when the surface energy of the gas-liquid interface is relatively small, the liquid will spontaneously spread in a thin film over the entire inner surface of the capillary tube. Now the capillary rise of the liquid in the tube is simply a matter of the pressure which is generated by a curved film of water that is tangential to the walls of the tube. In this case the radius of the curved surface is the radius of the tube wall and the pressure is given by Laplace’s equation for a spherical surface 42).

The Independent Domain Model of Hysteresis

If you confine yourself to the conceptual model presented in Figure (2.10), you will always be on either the main wetting or main draining curve, but Figure (2.8) shows that if one stops draining at some point short of θr, θ{Ψ} will depart from both the main wetting and main draining characteristic curves. This can be understood and accounted for through a simple extension of Figure (2.10). We noted that the key features of filling and draining of a single pore had to do with that pore’s neck and body radii, which we will refer to as rn and rb respectively. With this in mind, consider categorizing a medium by neck and body pore radius. As before, we will immediately translate these radii into equivalent capillary pressures via equation 26. So the idea is that we can characterize any medium as the aggregation of independent pores, each of which has a characteristic filling pressure (hf) dictated by the body radius and emptying pressure (he) controlled by the neck radius. [Although we appeal to equation 26 as a motivation, we can simply measure the filling and emptying pressures, so this is really just a crutch for visualizing the system.] Finally, we assume each pore is hydraulically connected to the bulk medium so that if a pressure is established at one of the medium’s boundaries, all the pores will experience that pressure. Each pore responds to changes in pressure independently (without regard to the state of the surrounding pores). This perspective was brought to vadose hydrology by A. Poulovassilis (1962).

The essence of this assumption is that the filling and draining of each pore is determined strictly by that pore’s geometry, regardless of the connection of that pore to surrounding pores (thus the term “independent”). This assumption is certainly fine for the wet end of the characteristic curve, when water fills most pores, but when the medium is quite dry it is possible the pores will become isolated hydraulically by being surrounded by empty pores. In this case the independent domain assumption breaks down for the isolated pore since it cannot drain regardless of the level of pressure exerted on the boundaries of the sample. We may expect then that the independent domain assumption will not lead to good predictions of residual moisture content, which must be carefully defined based on the process of interest, as discussed in section 2.2.

Figure 2.15
Figure 2.15. Cartoon of the characteristic radii for a pore and that pore’s characteristic filling and emptying pressures.

Pores will have a range of magnitudes for he and hf which can be described by a joint probability density function, as plotted in Figure 2.16. Note that the distribution functions for he and hf are not independent, since pores with very small necks are more likely to have similarly small bodies. The total volume under the curve in Figure 2.16 is 1, corresponding to the fact that all pores will have some combination of the two characteristic radii, which is stated mathematically in equation 68.

(68)   \begin{equation*} \int^{\infty}_{0} \int^{\infty}_{0} f \left(h_{e}, h_{f} \right) dh_{e} dh_{f} = 1 \end{equation*}

Pore necks are, by definition, smaller than pore bodies, so the absolute value of he is necessarily larger than hf for a given pore. So we may draw a line at 45° on the he-hf plane and know that all pores will lay below that line, as indicated in Figures 2.16 and 2.17.

Figure 2.16
Figure 2.16. Illustration of the probability density function moisture content as a function of filling and emptying pressures.
Figure 2.17
Figure 2.17. Identification of key features of the joint density function employed in the independent domain model of hysterisis given in Figure 2.16. The radius of the pore body is rb and that of the neck radius is rn.

To understand the density function plot more intuitively, it is useful to consider the physical interpretation of several of the points indicated in Figure 2.17. You see that the end-points of the pressure scales are defined by the largest and smallest bodies and necks. The appeal to Laplace’s equation here is primarily heuristic: the exact radii identified are not what you would find under a microscope, but provide a useful conceptual model. If we consider some horizontal line with a value h1 going from the right boundary (hemax) to the 45° line, this delineates all of the pores which fill at pressure h1. Since the filling pressure is determined by the body radius, we can consider these to be the pores with some body radius rb such that hf=2σ/r1. Similarly, the vertical line shown in Figure 2.17 at pressure he defines all the pores which empty at pressure he. Realizing that the pore neck diameter determined the draining pressure, we can again associate all of these pore as having neck diameters of rn such that he=2σ/rn. The intersection of these two lines corresponds to the pores which fill at hf and empty at he, and the height of f at this point indicates the fraction of pores which have this particular set of characteristic radii.

With this understanding of how to organize filling and emptying of pores, it is straightforward to obtain the various characteristic curves (once you have determined the joint density function). Applying this method boils down to figuring out the range over which to integrate the density function. As an example, consider the series of characteristic curves shown in Figure 2.8 within the context of the independent domain model. Starting from point 1 where the medium is dry, all the pores are empty. So now we will go from h = -∞ to h = 0. Through this sequence of pressures we want to add up all the pores which fill regardless of their emptying pressure. Thus we sweep downward through hf adding up all the pores with emptying pressures between h to -∞, integrating the density function as we go. This first wetting process is illustrated in Figures 2.18a and 2.18b.

Having saturated the sample, to obtain the main draining curve we follow the same procedure. This time the pressure starts at 0 and becomes more negative. Since all the pores are initially filled, we need to integrate along vertical lines all the pores which empty at a given pressure regardless of the pressure at which they filled (Figure 2.18c).

Figure 2.18
Figure 2.18. Mapping of the water filled pores through the sequence of wetting and draining characteristic curves presented in Figure 2.8.

Things get interesting when we stop prior to emptying all the pores, and turn around and start to re-fill the medium. Pores which are already filled with water cannot be refilled, which is the only feature which has changed from the main filling procedure. We still sweep down with a horizontal line, but now we only add the pores which fill between pressures h1 and hf (Figure 2.18d). Figure 2.18e continues this sequence, now following the draining process. So once the joint density function is known, it is simply a matter of determining which area of the domain to integrate over to determine the moisture content for any sequence of pressure changes. This procedure is quite amenable to numerical implementation.

This procedure can be stated mathematically as integrating over the domain of filled pores. We set up a few of the integrals which would be evaluated in Figure 2.8 as examples. In the main wetting process (Figure 2.18a), we simply carry out the integral

(69)   \begin{equation*} \mathrm{\uptheta\left(h \right)} = \int^{h}_{-\infty} \int^{h_{f}}_{-\infty} \mathcal{f} \left(h_{e}, h_{f} \right) dh_{e} dh_{f} \end{equation*}

Supposing we take the system to saturation, and now start draining the medium to obtain the main draining curve. Having filled all the pores to start with, we must calculate how much has been drained, and subtract that from 1:

(70)   \begin{equation*} \mathrm{\uptheta\left(h\right) = 1 - \int^{h}_{0} \int^{h_{e}}_{0} f (\left(h_{e}, h_{f} \right) dh_{f} dh_{e}} \end{equation*}

Suppose now that we stop draining at a pressure h1 and start refilling (Figure 2.18d). For convenience, let’s define the moisture content at the completion of the draining as

(71)   \begin{equation*} \mathrm{\uptheta\left(h_{1} \right) = \int^{h_{1}}_{-\infty} \int^{h_{f}}_{-\infty} \mathcal{f} \left(h_{e}, h_{f} \right) dh_{e} dh_{f}} \end{equation*}

In refilling, only the pores which were empty should be counted. So the integral will be restricted to the range of emptying pressures from hf to h6:

(72)   \begin{equation*} \mathrm{\uptheta(h) = \uptheta(h_{1}) + \int^{h}_{h_{1}} \int^{h_{1}}_{h_{f}} f(h_{e} h_{f}) dh_{e} dh_{f}} \end{equation*}

The only trick here is to keep careful track of the range of integration, with some bookkeeping: you have to keep track of the turning points of the soil’s history. By turning points, we mean the occasions when the soil either stopped wetting and began to dry, or vise-versa (i.e., there was a change in sign in the time derivative of the moisture content). A notation for recording the turning points was presented by Enderby (1955), and adopted by Everett (1955) and Mualem (1974)

(73)   \begin{equation*} \uptheta\left\{\mathrm{h} \right\} = \uptheta\left\{ \begin{array}{cc} \mathrm{h_{1}} & \mathrm{h} \mathrm{h_{0}} & \mathrm{h_{2}} \end{array} \right\} \end{equation*}

Figure 2.19. Experimental primary drying scanning curves (solid lines) and theoretical points obtained using the independent domain model (small circles). Also included is the joint density function for this medium. Reproduced with permission from Poulovassilis (1962).

In this notation, subscripts denote the order of pressures, and relative position indicates whether the transition was wetting or drying. In the case shown, the medium wetted from h0 to h1, dried to h2, and then re-wetted to the present pressure h.

In 73 we have written the relationship between pressure and moisture content using curly brackets {}. Later in this text, and in general in the literature, this relationship is written as a function with regular parentheses (), but this really is inappropriate, as θ is not a function of h but is a functional of h. That is to say, there is not a one-to-one mapping between θ and h without consideration of the antecedent conditions. Of course we can relate θ and h from a known initial state and through a known sequence of either θ or h as stated in 73, without requiring any knowledge of the other independent variables (most notably time). The ability of this model to reproduce hysteretic scanning curves is demonstrated in Figure 2.19 which is reproduced from Poulovassilis (1962).

So we have set up the conceptual framework for this method of modeling hysteresis, but how do we obtain the joint density function practically? You can just carry out a terrific number of experiments where you map out the entire domain of possible filling and draining pressures to obtain f(he,hf) by brute force. While this would appear to be impractical, with numerically controlled valves and scales, this is not infeasible using an automated pressure cell. However, in 1973 Mualem introduced a simplification of this model which greatly eased the task of determining the density function.

Mualem (1973) noted that the joint density function f(he,hf) could be well approximated by the product of two univariate density functions

(74)   \begin{equation*} \mathrm{f \left(h_{e} h_{f} \right) \approx g \left(h_{e} \right) l \left(h_{f} \right)} \end{equation*}

where g() and l() are probability density functions that depend only on he and hf, respectively. This constrains, for instance, the filling pressure distributions to be the same, up to a constant multiplier, along any draining pressure line, as illustrated in Figure 2.20. Using Mualem’s formulation, we only need the main filling and emptying curves to obtain g() and l(), and from these results all the hysteretic results can be derived. Mualem showed that this formulation maintained its ability to reproduce scanning curves (Figure 2.21).

Parlange (1976) went on to show that an even more constraining similarity model based on data from the main draining curve alone is sufficient to reproduce the full family of scanning curves. This was demonstrated to be capable of fitting experimental data by Hogarth et al. in 1988, as shown in Figure 2.22, reproduced from their work.

We have gone into a fair bit of detail in how to handle hysterisis, but have not mentioned when and why it must be included when describing vadose zone processes. First off, to understand the general character of many processes, dealing specifically with hysterisis is not critical. For instance in irrigated agriculture, as long as you know the moisture content at field capacity and the permanent wilting point, you may not worry about the exact pressure the soil achieves at intermediate water contents. On the other hand, hysterisis has been shown to be central to the nature of fingered flow processes (see Chapter 3), as well as the movement of NAPLs (e.g., Lenhard and Parker, 1987; Lenhard et al., 1991). Beyond these very specific settings the literature does not provide general guidance on when and where hysterisis is most important. Beyond the examples given, we really don’t know when hysterisis will play a central role.

If hysterisis has been shown to be important to several important processes, and in general we can’t be sure hysterisis will not be influential, why haven’t more people included consideration of hysterisis in the description of vadose zone processes? We think this is simply a matter of practical constraints. None of the analytical solutions to flow, and precious few of the numerical codes for modeling flow in the vadose zone consider the possibility of hysterisis. Until we have convenient numerical codes which include hysteretic effects, the possible effects of hysterisis will continue to be ignored in quantitative analyses. Of course nature doesn’t care if hysterisis is hard to model or easy to model: it happens! Failing access to a useful model which includes hysterisis, a responsible approach is to qualitatively assess the potential impact of hysterisis on a problem. Given the rapid increase in availability and sophistication of numerical models, we believe that hysterisis will be included in many popular models in the coming decade, and is in fact available in research models at this time (White et al., 1995).

Figure 2.20
Figure 2.20. Illustration of the implication of the similarity assumption introduced by Mualem (1973) for the independent domain emptying and filling density function.

Geometric Relations of Mualem Similarity:

Emptying pressure: \displaystyle{\mathrm{\frac{x} {y} = \frac{g (h_{4})} {g(h_{3})} = \upalpha}} and \displaystyle{\mathrm{f(h_{4}, h ) = \upalpha f(h_{3}, h)}} for any h

Filling Pressure: \displaystyle{\mathrm{\frac{z}{y}=\frac{l(h_{1})}{l(h_{2})}=\upbeta}} and \displaystyle{\mathrm{f(h,h_{1})=\upbeta f(h,h_{2})}} for any h

Figure 2.21. Mualem’s demonstration of the ability of the independent domain model with the similarity assumption to reproduce experimental scanning curves found by Poulovassilis (1970) (Reproduced from Mualem, 1974, with permission).

Figure 2.22. Primary draining curves as found experimentally and reproduced using the Parlange similarity model form of the independent domain model (Reproduced from Hogarth et al., 1988).

[Insert “Measurement of water potential” here]

Measurement of Water Content

The ability to accurately measure water content in time and space is central to the study of vadose zone hydrology. Remarkably, the ability to measure this critical parameter continues to present many challenges, with significant developments continuing to come forward.

The benchmark for water content measurements is the direct determination by oven drying (Gardner, 1986). In this method a sample of soil is gathered from the field and weighed. It is then placed in a 100 – 110℃ oven until it’s weight has stabilized. The mass basis water content may then be computed using 8. Though simple to describe, even this method is subject to significant error and interpretation. For instance, if stones are avoided, either as an unplanned feature of the method of collection or by design, the value of water content obtained will be systematically higher than the true field value. Furthermore, much of the water bound to clay remains in place at these temperatures. Thus you are not obtaining the total water content, but that which can be liberated under these specific conditions. Further discussion of these issues is presented by Gardner (1986). Perhaps the most problematic aspects of this method are that it is destructive and cannot be completed in the field in a timely manner which renders the method inappropriate for many applications.

There are a variety of field methods, all of which use indirect measurements of bulk soil physical properties to infer a moisture content. We will discuss three methods that illustrate the key physical properties that can be related to water content. These are the neutron thermalization method (neutron probe), the gamma ray attenuation method (gamma probe) and the dielectric method (time domain reflectometry).

Neutron Probe

The neutron probe has long been the most accurate method of measuring water content in the field. The method relies on the thermalization, or slowing, of neutrons. Relatively high energy neutrons are directed into the medium from a radioactive source, with a detector provided to observe the resulting number of slow neutrons that are back-scattered. These neutrons then collide with the various atoms in the soil. The loss of energy which occurs upon each collision is far greater when the collision is with an object of mass similar to the neutron. For instance, in a single collision the loss of energy from a neutron to hydrogen is about 4,000 time that of a collision with oxygen (Weinberg and Wigner, 1958). So the presence of hydrogen is closely correlated with the production of thermalized neutrons.

Since water is the largely hydrogen, the rate of thermalized neutron production can be calibrated to the surrounding volumetric water content. Site calibration is required to account for the local effects of the mineral and organic components and structurally bound water in the soil which will also contribute to neutron thermalization. The calibration equation is typically written in the form

(75)   \begin{equation*} \mathrm{\uptheta = a + b(I/I_{0})} \end{equation*}

where θ is the volumetric water content, Io is the incident neutron flux, I is the thermalized neutron count, and a and b are calibration parameters that depend upon the neutron source strength, soil properties, and detector characteristics. Since the soil characteristics will vary significantly in the features which will effect this calibration (e.g., organic matter), proper calibration must be carried out for each measurement location. This presents some logistical difficulties. First, the measurement must be taken at least two moisture contents to be able to obtain the two parameters. Therefore one must either come to the field twice, under differing moisture regimes, or one must carry out the calibration in the lab. Since the calibration requires a volume of soil at least as large as the instrument’s sphere of influence (about 15 cm in wet soil and 75 cm in dry soil), and in light of the sensitivity to bulk density, it is impractical to obtain sufficient samples to calibrate the probe using disturbed laboratory samples. This leaves only in situ calibration as being a viable option. This calibration is typically carried using the samples obtained from the original probe installation, then returning to the field when a second measurement can be made when the soil has quite different moisture content. In fact, while two points are sufficient in theory, a reliable calibration generally requires many more measurements, typically more like ten samples.

A great deal of work has been carried out to obtain estimates of the error associate with neutron probe measurements, calibration of measurements made close to the soil surface, and other logistical issues. A good starting place to learn about these practical concerns is Gardner (1986). Use of the neutron probe method has diminished in the past decade due to the logistical issues involved with handling and transporting radioactive sources. Yet due to the precision which may be obtained, the large sphere of influence of measurements, and low cost of operation, this method remains a mainstay of vadose moisture measurement.

Gamma Attenuation

The gamma attenuation method of liquid content measurement relies on the adsorption of gamma rays. The measurement is made by locating the sample under study between a gamma source and detector. The attenuation of the beam follows Beer’s law

(76)   \begin{equation*} \mathrm{I_{\uplambda} / I_{0\uplambda} = \exp(\sum \upalpha_{i\uplambda} x_{i}) = C \exp(\upalpha_{w \uplambda} x_{w})} \end{equation*}

where the subscript λ denotes the gamma frequency of interest, I is the detected gamma intensity without sample, Iλ is the detected gamma intensity after passing through the sample, α is the attenuation coefficient for material i, and xi is the path length of gamma travel in material i. If all of the features of the sample are constant except for liquid content, the attenuation of non liquid components can be identified with a single multiplicative constant C. For the fluid content, xi = Tθi where T is the beam path length in the sample.

Each gamma energy will be attenuated at differing rates for a given liquid. By using multiple frequencies of gamma radiation it is possible to track the content of multiple variable constituents. For this reason, the dual frequency gamma method is the most common method of documenting the joint movement of water/oil systems (e.g., Schiegg, 1990), and may be used to document swelling effects (Gardner et al., 1972).

Due to the need to place the sample between source and detector, the gamma method is not often convenient to apply in the field, but is widely used in laboratory environments. Furthermore, given the use of intense radioactive sources, the gamma method is limited to highly controlled settings. For an excellent overview of this method, see Gardner (1986).

Time Domain Reflectometry

The measurement of water content using time domain reflectometry (TDR), although long standing (e.g., Fellner-Feldegg, 1969), has rapidly grown in popularity with the advent of affordable TDR systems. The TDR method can be configured with multiplexing to allow nearly continuous automated measurement at multiple locations with a precision sufficient for most applications. Furthermore, TDR requires less calibration than the neutron probe method, and eliminates the use of radioactive sources. Finally, TDR can be used to measure bulk electrical conductivity in parallel with water content, which can be very useful in tracer and salinity studies. Disadvantages of the method are a relatively small sphere of influence (<10 cm radius), and the high cost of probes (on the order of $100 per measurement location). Due to the significant benefits of TDR and developmental nature of this technology, we will provide a slightly more in-depth discussion of this method.

TDR depends upon the measurement of travel time of an electromagnetic wave down a wave guide which is either inserted or laid on the soil surface (Selker et al., 1993; Maheshwarla et al., 1995). Travel time increases with the square root of increase in the real component of the dielectric surrounding the probe. Since the mineral portion of porous media typically have constant dielectric (between 2-5) at frequencies > 50 MHz (Hoekstra and Delaney, 1974; Campbell, 1990), and liquid water in such media has a dielectric of 70-80 (temperature dependent), variation in signal travel time along the probes is attributed to variations in water content. Typical accuracy of TDR measurements better than +/-2% volumetric (Roth et al., 1990).

With the travel times of signals along TDR probes being in the order of nanoseconds, the precise measurement of this travel time has only been economically feasible in the last 15 years with the advent of inexpensive dedicated TDR oscilloscopes. The specific implementation of the TDR method has been intimately tied to the characteristics of the measurement tools. In particular, the Tektronix 1502 series TDR meters which were designed for field evaluation of signal line integrity have been the primary instrument employed in TDR moisture measurement. A key determinant of spatial resolution in TDR is the rise time of the device; approximately 2×10-10 s for the Tektronix 1502 series. This dictates the minimum probe length that can be used. The change in travel time along a probe, Dt, for a change in travel velocity from v1 to v2 is given by

(77)   \begin{equation*} \mathrm{D_{t} = 2(L/v_{1} - L/v_{2})} \end{equation*}

where L is the probe length, with the factor 2 coming from the fact that two-way travel time is observed. Solving for L, we may calculate the minimum probe length for a given ability to measure Dt, and a given tolerance in velocity change. For example, when considering a dry soil-water system with a dielectric of 10, one can ask what minimum probe length allows determination of a change 2 percent in volumetric moisture content. This would create a change in dielectric of about 1.5, with v = 3×108 m/s x K1/2, the minimum desirable probe length is approximately 12 cm. Measurement accuracy can be increased, or the probe length decreased, when using a faster TDR system. Kelly et al. (1995) determined that probes as short as 0.02 m could be used with a TDR system which had a 25 ps rise time. The physical limitations of the frequency dependence of the dielectric of water become significant at these time scales. The period of molecular resonance of water is of the order of 20 ps. The real dielectric of water drops dramatically for signals with periods in this range, while the complex component (which is an energy loss term) climbs (Merabet and Bose, 1988). The combined loss of propagation delay and energy for single components with periods < 40 ps obviates any advantage of operating at frequencies in excess of about 10 Ghz. In addition to considerations of length, issues of spacing, conductor diameter, and number of wave guides are all significant in determining the sampling volume and sensitivity of the device (Knight, 1992; Whalley, 1993).

An additional limitation experienced in the application of TDR has been signal loss due to bulk conductivity of the porous media, which essentially shorts out the signal sent down the TDR probe (Arcone, 1986). The loss of signal strength due to the electrical conductivity of the media can be exploited to provide estimates of soil solution salinity, as demonstrated by Dalton et al. (1984) and many others. Signal loss difficulties are also encountered using Ground Penetrating Radar (GPR) methods in fine-textured media, as will be described later in this chapter. The signal loss problem can in part be mitigated through the use of minimum length insulated probes (e.g., Kelly et al., 1995), but it requires individual calibration curves for each probe geometry. In addition, high conductivity conditions reduces signal energy due to eddy currents established within the media, thereby limiting the range of applicability of the insulated probes to solutions with concentration of smaller than one molar total salinity, about 25 times greater than can be accommodated by un-insulated probes.

TDR methods are well established in vadose zone monitoring, yet there are various areas where significant enhancements can be achieved. One recent advancement is the development of multilevel shorting-diode TDR probe (Cook et al., 1992), which allows the monitoring of a linear moisture profile with the installation of a single probe. Another opportunity for enhancement of TDR methods lies in the interpretation of wave forms. At present, most TDR signals are interpreted with an empirical relationship between moisture content and travel time as developed by Topp et al. (1980). More recently, authors have derived physically-based methods to interpret the signals (Herkelrath et al., 1991; Heimovaara, 1993; Whalley, 1993), thereby allowing a more clear distinction between dielectric effects on the signal by the solid and liquid phase of the soil matrix. At this time, explicit electromagnetic models for signal propagation are being developed which predict the exact shape of the reflected wave form. These models can be used in a forward modeling mode, but the ultimate objective is to estimate phase contributions by inverse modeling using the full TDR trace as the basis of estimation of moisture content and soil conductivity, which will make more complete use of the wave form obtained from a TDR system.

Measurement of moisture content using TDR will continue to grow in the research community, but rapid adoption in the commercial sector has been hampered by the relatively high cost of TDR systems and the requirement of careful signal interpretation. New products have reduced these concerns, with increasingly robust signal interpretation software and more affordable meters. The TDR method has yet to be packaged into a system which is broadly accessible to individuals in the commercial sector. Nevertheless, TDR is continuing to grow in popularity in the measurement of moisture content within the vadose zone. Competing methods of moisture content measurement through sensing of the soil dielectric are under development or marketed (e.g., Campbell Scientific). As shown by Hoekstra and Delaney (1974) and Campbell (1990), these measurements must be made at frequencies of greater than 50 MHz to obtain a calibration which is independent of soil type.

Measurement of Hydraulic Conductivity

Hydraulic conductivity is one of the most important parameters in characterizing the transport of fluids and solutes through the vadose zone, but presents significant practical difficulty in measurement. In vadose applications, knowledge of the saturated conductivity is often of limited use, as these conditions are rare in most locations. Yet reliable methods of measuring unsaturated conductivity are very few. One of the most common approaches to obtaining the unsaturated conductivity relationship is to measure both the saturated conductivity and the water retention curve, and estimate the unsaturated conductivity based on a conceptual model relating permeability to water content (van Genuchten, 1980). This process is facilitated through application of computer programs written expressly for this purpose, such as RETC (van Genuchten et al., 1991).

The only direct method of measuring unsaturated conductivity in the field is through a pressure/moisture content profile method. The most widely used approach is the “instantaneous profile method” (IPM). The IPM requires that an area of land be subjected to ponded infiltration until the entire profile of interest is saturated. The site is instrumented with devices to measure both the moisture content and pressure with depth and time. At t = 0 the ponding is ceased, and the site is covered with an impermeable barrier to avoid evaporative losses. Pressure and water content measurements are made in time, with measurements extending to periods of weeks and months if low permeability values are required. Using the water content measurements one can calculate the flux at various depths by way of the total time rate of change in water content in the profile above the height of interest. By observing the local pressure gradients, as measured with tensiometers, one then has values for q and dh/dz in Darcy’s law, and it is straightforward to solve for K. A very nice added benefit of this method is that the joint pressure/water content measurements provide an in situ determination of the water retention relationship. Although the IPM is a direct measurement, due to the large number of measurements that need to be made over long periods of time, it is difficult to carry out, and often provide results which are less than satisfactory (e.g., Selker, 1994).

The remaining methods of measurement rely on an analytical solution for the flow from some simple geometric boundary condition. The solution must predict the time rate of infiltration in terms of the desired parameters. The non-linearity in Richards equation generally leads to the adoption of a solution with very restrictive assumptions in order to allow development of a general solution. The two most widely used assumptions are (1) that the infiltration process may be approximated by a sharp wetting front followed by a region of nearly saturated media (Green and Ampt infiltration: see Chapter 3) and (2) that the hydraulic conductivity function may be approximated by the Gardner (1958) relationship

(78)   \begin{equation*} \mathrm{K = K_{s} \exp(\upalpha h)} \end{equation*}

where α is a parameter approximately equal to -1/hwe, and thus is large for coarse textured soils (on the order of 10 m-1 for sands) and small for fine textured soils (on the order of 0.1 m-1 for clays).

The best control over the boundary at which infiltration occurs is obtained by carrying out surface measurements. Ring infiltration methods cover a broad set of these tests. In general these method involve establishing a prescribed pressure over a circular area and observing the resulting infiltration. The scale of these measurements ranges from 100 cm2 (Ankeny et al., 1988) to 4 m2 (ASTM, 1990). An approach that make use of the Green and Ampt concepts was introduced by Bouwer, and is discussed in section Here we will introduce the ponded double ring method, and some newer tension infiltrometer approaches.

By far and a way the most common hydraulic test carried out for to estimate soil hydraulic properties is the double ring infiltrometer (Figure 2.40). The rings are inserted deep enough to preclude leakage from the outer ring, and to have the tops of the rings level with each other. A constant water level is quickly established in both rings to the same level, and the infiltration of water may be measured by watching the drop in water level on a floating ruler, or by using a Marriott bottle constant head water source.

The data may be analyzed using any of several infiltration models. Below we list the Horton Equation, the two-term Philip model, the Green and Ampt model and the Brutsaert model.

Horton Equation (1940): In this empirical infiltration model the rate of infiltration, i, is given by

(79)   \begin{equation*} \mathrm{i = i_{f} + (i_{0} - i_{f}) \exp(-\upbeta t)} \end{equation*}

where if is the infiltration rate after long time, io is the initial infiltration rate and β is and empirical soil parameter. Integrating this with time yields the cumulative infiltration

(80)   \begin{equation*} \mathrm{I = i_{f}t + \frac{i_{0} - i_{f}} {\upbeta} (1-\exp(-\upbeta t))} \end{equation*}

The model is widely used, but since the parameters obtained are not directly associated with physically measurable soil properties (e.g., conductivity), it is of limited use except for predicting infiltration.

The two term Philip model (1957) is a simplification of the infinite series Philip introduced for infiltration, and is widely used. The result given by

(81)   \begin{equation*} \mathrm{i = 0.5 S t^{1/2} + A} \end{equation*}

or for cumulative infiltration

(82)   \begin{equation*} \mathrm{I = S t^{1/2} + At} \end{equation*}

has the advantage of being closely related to measurable soil physical parameters (i.e., S is the soil sorptivity and A is typically taken to be an estimate of field saturated conductivity).

The Green and Ampt Model (1911) for constant head infiltration is derived in Chapter 3, and is also useful in interpreting infiltration results. The Green and Ampt result for constant head infiltration is given by

(83)   \begin{equation*} \mathrm{t = \frac{n}{K_{sat}} \left\{L - (h_{f} + d) \ln \left[\frac{h_{f} + d + L} {h_{f} + d}\right] \right\}} \end{equation*}

where L is the depth of wetting front penetration, n is the porosity of the soil, d is the depth of ponding, and hf is the water entry pressure for the soil. The cumulative infiltration is simply I = nL. To use this equation you must find the values of Ksat and hf which give the best fit to the data.

The Brutsaert Model. Brutsaert (1977) revisited the infinite series solution of Philip (1957) and developed a very convenient result that provides an excellent approximation to the full series solution for infiltration, and is written explicitly in the physical variables of interest. For the rate of infiltration Brutsaert (1977) found

(84)   \begin{equation*} \mathrm{i=K_{sat} + 0.5 S t^{-1/2} \left[1 + \upbeta \left(\frac{K_{sat} t^{1/2}} {S} \right)\right]^{-2}} \end{equation*}

where S is the sorptivity and 1/3<β<1 is a soil parameter related to the distribution of pore sizes. For soils with very wide pore size distributions you should take β close to 1, for most other soils β = 2/3 is quite good, while for uniform pore size soils, β = 1/3 is suggested. The cumulative infiltration using the Brutsaert model is simply

(85)   \begin{equation*} \mathrm{I = K_{sat} t + \frac{S^{2}} {\upbeta K_{sat}} \left\{1 - \left[1 + \upbeta\left(\frac{K_{sat} t^{1/2}} {S} \right)\right]^{-1} \right\}} \end{equation*}

Fitting this equation to infiltration data allows ready determination of both Ksat and S.

Figure 2.40
Figure 2.40. Typical double ring infiltrometer installation. Ring diameters employed vary widely.

At negative pressures, features of textural interfaces and fissures cause drastic reduction in unsaturated permeability, while at near-zero pressure macroporosity can greatly increase permeability. Thus, the ability to measure hydraulic properties of vadose materials in-situ at realistic negative pressures has been recognized as an area in need of attention.

Controlled infiltration under tension to determine unsaturated water retention and transport properties has been considered for decades, however, only in the past few years have experimental and analytical tools been improved to make the technique practical in the field. Clothier and White (1981) introduced the “sorptivity tube” which allowed for fixed pressure (<0) application of water to a soil using a Marriott bottle supply with a sintered glass porous interface. This design was tested and improved, culminating in the fabric interface device presented by Perroux and White (1988). This device provided ready adjustment of interface pressure, could be made with a wide range of contact areas, and was robust and simple enough to be used conveniently in the field. Developed in parallel, Ankeny et al. (1988) fitted the tension infiltrometer with two pressure transducers which allowed automated acquisition of both the pressure at the soil interface and the flux from the device. Further, this device was fitted with three adjustable bubble tubes that allowed for rapid evaluation of flux from a series of interface pressures. A key development in interpretation of the results of these devices was identification of the utility of the result of Wooding (1968) to disk infiltrometer (Ankeny et al., 1991). Wooding found that using a Gardner exponential conductivity function 78 the flux into a circular wetted area of radius r, at steady state, was

(86)   \begin{equation*} \mathrm{Q = \uppi r^{2} K_{s} \exp(\upalpha h)\left(1 + \frac{4} {\uppi \upalpha r}\right)} \end{equation*}

where Q (measured) is the flux into the soil over a device of radius r (measured), h (measured) is the pressure at which the infiltration boundary, Ks (unknown) is the saturated hydraulic conductivity, and α (unknown) is as defined in ??. By carrying out infiltration tests at two or more pressures, or with two or more radius disks, one may solve for both of the unknowns in 86.

One significant shortcoming of Wooding’s approach is the notion of a well-defined steady-state condition, and the time required to achieve this condition. As shown by Smettem et al. (1995), these effects may take hours to become small, and the outflow might appear to be steady far before the true steady state has been achieved. Haverkamp et al. (1994) present an approximate solution for the time-dependent rate of infiltration for one and three dimensional cases as

(87)   \begin{equation*} \mathrm{I_{1d} = S_{0} t^{1/2} + \left(K_{n} + \frac{2 - \upbeta} {3} (K_{0} - K_{n}) \right)t} \end{equation*}

(88)   \begin{equation*} \mathrm{I_{3d} = S_{0} t^{1/2} + \left(K_{n} + \frac{\upgamma S_{0}}{(\uptheta_{0} - \uptheta_{n})r} + \frac{2 - \upbeta} {3} (K_{0} - K_{n}) \right)T} \end{equation*}

where So is the sorptivity, β is a shape parameter (0<β<1), γ is a parameter with value between 0.6 and 0.8 which Smettem et al. (1994) found to be well approximated by 0.75, Ko & Kn, and θo & θn are the conductivity and moisture content at the antecedent and instrument pressures respectively. As shown by Smettem et al. (1995), a double ring tension infiltrometer may be made in which the outflow of the inner ring may be assumed to be one-dimensional, and that of the combined inner and outer rings three-dimensional.

Smettem et al. (1995) show that a good estimate of sorptivity may be obtained using 88 by plotting I3d/t1/2 against time. We see that 88 then reduces to

(89)   \begin{equation*} \mathrm{I_{3d} / t^{1/2} = S_{0} + At} \end{equation*}

for some constant A. By fitting a straight line to 89, the t=0 value of this relationship provides an estimate of So. Smettem et al. found this value to be consistent with values obtained from one-dimensional infiltration tests. This result may be used in conjunction with 86 or 88 for Ks.

These measurements are typically made at a series of tensions and often with two sizes of circular sources. The measurements require access to a horizontal soil surface, thereby limiting application of this technique to the first few meters of the vadose zone, although a down-hole tension infiltrometer designs has also been developed (Sisson and Honeycutt, 1994).

An alternative method of interpreting tension infiltrometer results is to carry out the estimation of sorptivity and hydraulic conductivity functions by inverting a numerical simulation of the process (Simunek et al., 1995 and 1996). Since this method makes full use of the early time data, it holds the promise of significantly increasing the utilization of site-specific hydraulic properties by reducing the experimental time required to obtain reliable values. The accuracy and practicality of these methods has also been explored recently, comparing a number of methods both for field and laboratory application (e.g., Logsdon and Jaynes, 1993; Gupta et al., 1993; Stolte et al., 1994; Cook and Broeren, 1994). It is clear that each of the present implementations of the tension infiltrometer method have difficulty under some situations, but the outlook for improved methods of measurement appears quite bright in the coming years using this technique.

To determine the hydraulic properties of vadose materials at depth, down-hole methods are essential. This is an area of active development, with some very useful techniques available. The method most widely employed is to auger a hole of constant radius to the depth at which the measurement is required. A constant head of water is then established in the hole, while the rate of infiltration is observed. Several tools have been developed for this test, and generally consist of a Marriott bottle source of water with an extension tube to reach the required depth. Such tools can easily be made without access to specialized equipment, as shown in Figure 2.41.

Figure 2.41
Figure 2.41Figure 2.41. Simple down-hole constant head permeameter. The graduated cylinder is filled prior to the test, then the device is lowered into the soil until the support plate sits on the soil surface. The test is commenced by opening the shut off valve, at which time the bubble tube establishes a zero pressure level at its point of entry to the drop pipe.

The analysis of constant head well permeameter can be carried out using a number of formulations. The original work directed to analysis of such data was presented by Glover (Zangar, 1953). As presented by Elrick and Reynolds (1993), Glover found that at steady state the rate of infiltration from such a permeameter could be approximated by

(90)   \begin{equation*} \mathrm{Q = 2 \uppi H^{2} K_{fs} / C} \end{equation*}

where Kfs is the “field saturated” hydraulic conductivity and C is a dimensionless shape coefficient reflecting the geometry of the wetted area. This was estimated by Glover to be

(91)   \begin{equation*} \mathrm{C = \sinh^{-1} (H/a) - [(a/H)^{2} + 1]^{1/2} + a / H} \end{equation*}

This solution assumes complete saturation in the wetted region, and does not include gravity driven flow (Elrick and Reynolds, 1993). For ratios of H/a < 2, this solution has been shown to be reasonably accurate, while for H/a > 10, it under estimates the correct value of C by about 50%, leading to an error in Kfs of about 40%. To account for the effects of gravity and unsaturated flow effects, Elrick et al. (1989) provided the solution

(92)   \begin{equation*} \mathrm{K_{fs} = \frac{C Q}{2\uppi H^{2} + \uppi a^{2}C + 2\uppi H / \upalpha^{ast}}} \end{equation*}

Elrick and Reynolds (1992a and 1992b) suggest that α* may be estimated based on field conditions using the values given in Table 2.5. Alternately, tests using multiple depths of ponding may be carried out, providing for multiple equations from which the values of Kfs and α* may be estimated using regression analysis. Difficulties with this approach have been reported (e.g. Elrick and Reynolds, 1992b) due to heterogeneity in soils which will yield inconsistent values of the parameters depending upon the depth of ponding. This can lead to non-physical values of estimated parameters (e.g., negative conductivity), and thus should only be employed in soils without significant vertical variability in properties.

Table 2.5 Suggested values for α* as provided by Elrick and Reynolds (1992).

α* Comments
Essentially the Glover Solution. Probably applicable to gravels.
36 Coarse sands and highly structured soils.
12 Most structured soils and medium and fine sands.
4 Unstructured fine-textured soils.
1 Compacted clays (e.g., clay liners).
0 The Gardner solution (Reynolds and Elrick, 1985). Pressure and gravity contributions negligible.

Review Questions

  1. What are the units of surface tension and why?
  2. What is the definition of the reduced radius of curvature and what is the pressure drop across a curved interface?
  3. If the local atmosphere in a soil is, at equilibrium, below the saturated vapor pressure for water, is there any liquid water for plants?
  4. What is a Haines Jump?
  5. What is the raindrop effect?
  6. What is hysteresis?
  7. What are the two basic assumptions of the independent domain model for hysteresis?
  8. Is Darcy’s law linear? in which terms?
  9. When does Darcy’s Law break down? What are the Reynolds number criteria?
  10. Is the unsaturated version of Darcy’s law hysteretic in moisture content? Pressure? Why?
  11. What are the two relationships that go into the derivation of Richards equation?
  12. Is Richards equation linear in any terms? (Consider K, θ, D, h, x, z)
  13. What is the point of using the D parameterization instead of K in Richards equation?
  14. Why have so many people developed functional relationships for D and K?
  15. What are the basic assumptions of Miller Similarity?
  16. List the dependence of the following terms on the characteristic grain size of the media (don’t worry about the effects of viscosity and surface tension): K, S, h.
  17. What is a NAPL, an LNAPL, and a DNAPL?
  18. Explain why the observation well depth of a NAPL is a poor indicator of true “free NAPL” content at a site.
  19. What is a tensiometer, and how does it work?

2.12 Problems

  1. The raindrop effect illustrates hysteresis in fluid/solid contact angle.
  2. A Haines jump is an Olympic sporting event developed by Sir Walter Haines in the late 1920’s.
  3. For a soil at equilibrium, if the soil gas has a relative humidity of 0.5, the moisture content is approximately 0.5 θs.
  4. Darcy’s law breaks down for high flow when the flow regime becomes non-laminar at a Reynolds number of about 100. Below this value flux and gradient in potential are linearly related.
  5. Darcy’s law breaks down for Reynolds number greater than 1 due to turbulence at the pore scale.
  6. Since conductivity is hysteretic in moisture content, even if you know the local gradient in pressure, it is not possible to be sure of the local flux at a given moisture content without knowledge of the history of the wetting at this point.
  7. Under steady state conditions, Richards equation reduces to Darcy’s law for unsaturated flow since conservation of mass is satisfied.
  8. The units of surface tension are energy per area.
  9. Poiseuille’s law tells us that conductivity follows a linear relationship with average pore size.
  10. The Green and Ampt infiltration model would be applicable to rapid infiltration processes for a falling head situation.
  11. A soils sample taken from the field with a volume of 42 cm3 weighed 29 grams, when saturated weighed 32 grams, and when dried weighed 21 grams. Calculate the soil’s moisture content (both volumetrically and on a mass basis), saturated moisture content, porosity, bulk density, particle density, and void ratio.
  12. For a fluid with surface tension of 2 N/m, and with curvature defined by the two radii -0.02m and -0.01m, what is the pressure in the liquid if the gas pressure is 2N/m3? Express your result in units of N/m2.
  13. Assume that we have a medium made up of pores identical to the one shown below. Assume that the total porosity of the medium is 28%. Plot the wetting and draining curve for the medium in the axis shown below. Assume a surface tension σ and a contact angle of γ. (You may assume that L >>R).13a13b
  14. A soil is found to respond as if it were made up of a capillary tube bundle as shown on the right. For each capillary of radius R there are two tubes of radius R/2, and four tubes of radius R/4. If the saturated conductivity of the medium is Ks, calculate the unsaturated conductivity and pressure potential of the medium for the three possible values of θ (assume that there is no hysteresis and write the range of Ψ in units of 2σcosγ/R, θ in units of θ/θs and K(θ) in terms of K/Ks). Put the results in the table given.
    Moisture Content

    (units of θ/θs)


    (units of K/Ks)

    Range of h

    (units of 2σcosγ/R)

  15. Solve for and plot the values of k(h) for h =-9.2, -26.7,- 49.4, -69.1, and estimate k(-90 cm) using the data obtained in a steady-state, homogeneous, vertical column experiment shown below.
    Pressure (cm) Elevation (cm)
    0.0 0
    -9.2 10
    -18.2 20
    -26.7 30
    -34.9 40
    -42.4 50
    -49.4 60
    -55.5 70
    -60.9 80
    -65.4 90
    -69.1 100

    K(-9.2) = _______________

    K(-26.7) = _______________

    K(-49.4) = _______________

    K(-69.1) = _______________

    K(-90) = _______________

  16. Plot the characteristic scanning curves for θ(h) for the normalized f(he,hf) data and the following wetting and drying history.

        \begin{equation*} \mathrm{\uptheta(h) = \uptheta} \left\{\begin{align}&-40& &-30\\0& &-10&\end{align}\right\} \end{equation*}


  17. Relative Humidity in the soil is useful for calculating pressure head at high suctions using the psychrometric equation. In semi-arid climates, soil water pressures may exceed -50 bars (about -750 psi, or -500 m of head!). Make a graph of relative humidity vs. pressure head over the range of 0 to -50 bars.
  18. Porous filters can be obtained with average pore size of 0.2 microns in radius. [Of course, the largest pores are those that dictate the air entry pressure of the plates.] Such porous plates can be used to de-water soil at high pressures to obtain the pressure-moisture content relationship. Assuming that the largest continuously connected pores through the plate are 2.5 times this size, and that the plate is initially completely saturated with water, at what pressure will the plate let air pass through? Express your answer in units of N/m2, psi, Bars, atmospheres, and cm of water.18
    1. Calculate the capillary pressure in kPa for the water shown between the two spheres (State your assumptions).
    2. What is the force between the spheres? [Think about balancing the forces of the surface tension, the contact force between the spheres and the pressure of the water. It is easiest to calculate about the plane of symmetry.]
    3. Why does wet sand stick to your hands? What about making sand castles with dry sand vs. wet sand?
      The scale of the drawing is 1 cm = 0.1 mm.
  19. Another way to obtain the Laplace equation is to balance the forces between surface tension and internal pressure for a sphere floating in space. Carry out this calculation for the force balance across the middle of a sphere of water in air. Calculate the pressure in a 0.1 mm sphere of water.
  20. Characteristic draining curves(adapted from similar problem by D.B. Stephens).

    Yechezkel Mualem has done a great service to folks studying the flow through the vadose zone with the publication of the text A Catalogue of the Hydraulic Properties of Unsaturated Soils in 1976. This catalogue collects the characteristic draining and wetting curves, bulk density, particle size distribution, and conductivity. You will find two typical entries from this catalogue attached; one for “Gilat Loam”, and the other for “Sinai Sand.” Consider a soil profile as shown below, which might have resulted from a glacial deposition (In solving this problem, ignore the effects of hysteresis and express pressure in the units cm of water).

    1. Plot the characteristic draining curves for the two soils using Mualem’s data.
    2. Plot the pressure versus elevation for the profile for the transects labeled A-A’ and B-B’ assuming that the profile is at steady state with the water table at 150 cm below the soil surface, as shown.
    3. Plot the moisture content versus elevation for the transects labeled A-A’ and B-B’.
    4. Suppose that the soil started by being saturated to the soil surface, and then was drained to have the water table at a depth of 150 cm, as shown. How much water (in cm depth) drained from the profile?
  21. Hysteresis Curves using Independent Domain Model. (after problem by Dr. W. Brutsaert, Cornell University).

    Using the data for the pore space density for Adeline Dune Sand (Talsma, 1970) given below, draw the pressure vs. moisture content for the following pressure series:

    (93)   \begin{equation*} \mathrm{q(y) = q \left\{} \begin{align} &0& &-30&\\ -100& &-70& & -50 \end{align} \right\} \end{equation*}

    water content units: 10-4(cm3 water/cm3 soil)/(10 cm water)2
  22. The diameter of a soil cylinder is 0.12 m and its height is 0.28 m.

    (i) What was the gravimetric water content (θm) if the wet weight was 5.5 kg and the dry weight was 4.25 kg? (ii) What was the soil bulk density ρb (assume particle density of 2650 kg/m3); and (iii) what was the volumetric water content θv?

  23. What were the porosity (n), the void ratio (e), and the degree of saturation (S) of the soil sample in problem 5?
  24. How many 200-liter barrels of water would be stored in the top 1.2 m of 1 hectare of land if the average θv was 0.25? What would be the volume of air-filled pores if the porosity (n) was 0.45? How deep would 80 mm of water wet a soil with ρb=1250 kg/m-3 and initial θv=0.14 m3/m3? Assume complete saturation immediately after wetting, and ρs=2650 kg/m3. What would the final depth of wetting be if the soil then drained internally to a field capacity water content of θv=0.26?
  25. Derive the following relations:

    (i) \mathrm{e = n/(1-n)}; (ii) \theta_{\mathrm{v}} = \mathrm{S n}; (iii) \mathrm{n_{a} = n (1-S)} where \mathrm{n_{a} = V_{a} / V_{t}} is the air-filled porosity.

    Example: derive \mathrm{n} = 1 - \rho_{\mathrm{b}} / \rho_{\mathrm{s}}

        \begin{align*} &\mathrm{n} = \frac{\rho_{\mathrm{s}} - \rho_{\mathrm{b}}}{\rho_{\mathrm{s}}} = \mathrm{\frac{\frac{M_{s}}{V_{s}} - \frac{M_{s}}{V_{t}}}{\frac{M_{s}}{V_{s}}} = 1 - \frac{V_{s}}{V_{t}} = \frac{V_{t} - V_{s}} {V_{t}}} \\ &\Rightarrow \quad \text{since } \mathrm{V_{t} - V_{s} = V_{f}} \text{ we obtain the definition of porosity as} \quad \mathrm{n = \frac{V_{f}}{V_{t}}} \end{align*}

  26. The following measurements were collected for the calibration of a neutron probe:
    Count Ratio – CR 1.7 1.5 1.2 1.0 0.8 0.6 0.4
    Water Content qm [kg/kg] 0.27 0.23 0.17 0.15 0.11 0.09 0.06

    The Bulk density was 1300 kg/m3. Determine the calibration equation for the neutron probe and plot the data with the calibration line.

  27. The distance to a TDR trace initial reflection (x1) was 5.01 m and to the final reflection (x2) 5.63 m. Given a probe length of 0.2 m and Vp=0.99; (i) What is the bulk dielectric constant εb=[(x2-x1)/(L Vp)]2, and (ii) Using Topp’s equation, determine the volumetric water content of the soil.
  28. For a soil with ρb=1300 kg/m3 what would be the dielectric constant (εb) at θv=0.22? (hint: use Eq.18 ).
  29. What would be the volume of drainage flowing into 1.5 m deep drains from a 2.15 hectare field with θvFC=0.22 m3/m3 and initial θv=0.18 m3/m3 following 150 mm of rain?
  30. Given a depth-dependent field capacity as: θvFC=0.26-0.03z, where z is the depth in m, and assuming a constant wilting point at θvWP=0.04 m3/m3, compute the volume (m3) of plant-available soil water in the top 1.5 m of a 25 hectare field.
  31. What is the gauge reading of a tensiometer if the matric potential ψm = -2500 mm and the gauge is 300 mm above the porous cup?
  32. What is hydraulic potential ψh for a situation where ψm = -280 mm and the reference elevation is 520 mm below the point of interest?
  33. The cups of tensiometers 1 and 2 are at a depth of 0.6 and 0.8 m, respectively. The gauges are 0.2 m above the soil surface. The gauge in tensiometer 1 indicates ψgauge = -0.9 m. (i) Draw the potential diagram (use weight basis or head) assuming static equilibrium conditions (use the soil surface as reference level). (ii) What is the gauge reading in tensiometer 2?
  34. If the EC of a soil solution was 21 dS/m, and ψm = -280 kPa, what was the relative humidity of the soil air at 20℃?
  35. What should be the effective diameter of the largest pore in a tensiometer cup used to measure pressure heads of -8.0 m?
  36. What is the height of rise in clean glass capillaries having diameters of 0.01, 0.05, 0.3, and 1.5 mm, when dipped vertically into a pool of pure water?
  37. What would be the effect on the height of capillary rise of problem 27 if the capillaries were tipped vertically into a pool of mercury, with the contact angle mercury-glass being 148°. Explain the differences.
  38. Use the SWC and van Genuchten parameters α=0.356 m-1, n=1.501, θr=0.039 m3/m3, and θs=0.418 m3/m3 to estimate plant available soil water between field capacity (ψm=-3 m) and wilting point (ψm=-150 m). What are the values of θvFC and θvWP?
  39. Fit the van Genuchten and Brooks and Corey parametric models for Soil Water Characteristics to the following experimental data.
    Matric Potential [-m] 0.01 1 2 3 5 8 61 110 323
    Water Content [m3/m3] 0.433 0.347 0.300 0.263 0.222 0.183 0.116 0.093 0.066
  40. What is the volume (m3) of plant-available soil water stored in the top 1 m of a 2 m2 area of Sarpy loam soil that is characterized by the following SWC VG model parameters: α=2.79 m-1; n=1.6; θs=0.4 m3/m3; and θr=0.032 m3/m3
  41. The relative humidity measured with a psychrometer at 25℃ is 0.45 (45%). What is the water potential yw of a soil sample in equilibrium with ambient air?
  42. A uniformly packed, U-shaped soil column has one end suspended in a container of water as shown in the figure at the right. The water level in the container has been kept constant and evaporation has been prevented from the open end for a period sufficient that the soil column has reached equilibrium. Thus there is no net flow in the column.

    Find: The hydraulic potential ψh and its components ψp, ψm, and >ψz for points A through F in the soil column. List your results.




atmospheric pumping


Beer’s law

Brutsaert Infiltration model

Buchner funnel


Bukner funnel

bulk density



capillary tube

characteristic curves

Characteristic Curves

chemical potential

conductivity, down hole measurement

constant head permeameter


contact angle


creeping flow

curved fluid surface


Darcian Length

Darcian velocity

Darcy’s law


dead end pore


dew point






double ring infiltrometer

down-hole tension infiltrometer


dynamic viscosity





field capacity

field saturated

Fokker-Planck equation

free product


fundamental differential form


gamma attenuation

gamma method, dual frequency

Gardner conductivity

gas movement

gas permeability


geometric similarity

Gibbs-Duhem relationship


gravitational potential

Green and Ampt

ground penetrating radar

gypsum block


Hagen Poiseuille

Hagen-Poiseuille Equation

Haines jump

Henry’s law

Horton Equation

hydraulic conductivity, measurement

Hydrostatics, 21

hysteresis, gypsum blocks





independent domain assumptions

inertial forces

infiltration models

ink bottle effect

instantaneous profile method

internal energy

intrinsic permeability



kinematic viscosity

Klinkenberg effect


Liquid Content


LNAPLs in wells



matric potential

mean free path length


Miller similarity




Navier Stokes Equation

negative pressure

neutron probe

neutron probe, errors

Newtonian fluid

non aqueous phase liquids

non-spreading NAPL

non-wetting fluid


oil spill


osmotic potential

oven drying


particle size distribution


permanent wilting point



pore size distribution


potential energy

pressure head

pressure potential

pressure potential, measurement


psychrometric equation



rain-drop effect

reduced radius

relative humidity

Representative Elementary Volume



Reynolds number

Richards equation




saturated liquid content

scaled conductivity

scaled pressure

scaled sorptivity

scaled time

scanning curve

similar media


similarity assumptions

soil diffusivity

sorptivity tube

specific flux

specific gravity

spreading NAPL


surface tension

Surface Tension

surface tensions



TDR, insulated probes

TDR, shorting diode



tension infiltrometer

tension infiltrometer, numerical simulation




time domain reflectometery

Time domain reflectometery


Total Bulk Density

total head


two-term Philip model




Van der Waals

vapor extraction

vapor pressure

Vapor Pressure

viscous dissipation

viscous forces

void fraction

Void Ratio

volatile organics


water content, measurement


Wooding solution


Youngs-Laplace equation



Brown and Montgomery, 1913, “Dehydration of Clays,” Bureau of Standards Technologic Paper no. 21, Washington, D.C.: https://nvlpubs.nist.gov/nistpubs/nbstechnologic/nbstechnologicpaperT21.pdf

Kirkham and Powers (1972)

Revil, A. and A. Jardani, The Self potential Method: theory and applications in environmental geosciences, Cambridge University Press, ISBN 978-1-107-01927-0, 2013.

Sposito (1989)

†† A few words about the terms in the K tensor.  So the idea here is to be able to represent the fact that the resistance to flow might be different per direction due to orientation of the particles at a scale below that of the REV.  So we will “take apart” the conductivity to explicitly represent the relationship between the gradient in potential in a particular direction, and the resulting flow that might happen in that same direction, or along a different axis.  Writing Darcy’s law

    \begin{align*} \mathrm{q} &=- \left[ \begin{aligned} \mathrm{{K_{xx}} {K_{xy}} {K_{xz}}} \\ \mathrm{{K_{yx}} {K_{yy}} {K_{yz}}} \\ \mathrm{{K_{zx}} {K_{zy}} {K_{zz}}} \end{aligned} \right] \left[ \mathrm{\frac{\partial h} {\partial x} \frac{\partial h} {\partial y} \frac{\partial h} {\partial z}} \right] \\ \\ &= -\left[ \mathrm{{K_{xx}} \frac{\partial h} {\partial x} + {K_{xy}} \frac{\partial h} {\partial y} + {K_{xz}} \frac{\partial h} {\partial z} ; {K_{yx}} \frac{\partial h} {\partial x} + {K_{yy}} \frac{\partial h} {\partial y} + {K_{yz}} \frac{\partial h} {\partial z} ; {K_{zx}} \frac{\partial h} {\partial x} + {K_{zy}} \frac{\partial h} {\partial y} + {K_{zz}} \frac{\partial h} {\partial z}} \right] \end{align*}

To get a feel for the terms in the K-tensor, look at the terms in the x-coordinate of flux.  We see that gradients in each of the x, y, and z directions can give rise to flux in the x direction.  In general we see that the term Kab is the proportionality constant relating gradients in potential in the b-direction to flux that results in the a-direction.  How could such a thing happen?  If the medium is aniostropic (i.e. permeability is a function of direction) then gradients that are not aligned with the directionality of permeability will give rise to flux not parallel with potential gradients.  An easy analogy is movement on a slide in a playground.  The direction of energy gradient due to gravity is vertically downward, but your resulting motion includes a lateral component, along the direction of the slide.  This really is similar to the motion of water in a clay with oriented particles:  the water will largely travel parallel to the sheet-structures of the clay rather than across the clay sheets because the resistance is so much lower, even if the gradient is not exactly in the same direction as the orientation of the clay sheets. If the coordinate system is aligned such that the coordinate system and primary directions of anisotropy are coincident the off diagonal terms in the permeability tensor will be zero (i.e., the Kab terms where a ≠ b), leaving permeabilities in the three coordinate directions.  If these are all equal, then the tensor collapses to a scalar.  So the only reason to take the permeability as a tensor is to capture the effects of anisotropy.

Some of you may ask what is this up-side-down triangle all about?  The ∇ symbol stands for the “dell” operator.  It is just mathematical short hand for calculating the three-directional derivative instead of the one directional derivative.  It is defined as

    \begin{equation*} \mathrm{\nabla = \frac{\partial} {\partial x}i, \frac{\partial} {\partial y}j, \frac{\partial} {\partial z}k} \end{equation*}

where i, j, and k are the unit vectors in the x, y, and z directions.  The result of “operating” on a function with ∇ is to obtain a new vector which is the slope of the function.   In practice, for any scalar valued function of x,y,z, ∇F tells us the magnitude and direction of the slope of F: ∇F points directly towards the steepest direction of uphill with a length proportional to the slope of the hill.   Later you will also see this symbol in the notation ∇·F.  The dot just tells us to take the dell and calculate the dot product of that and the function F (which needs to be a vector for this to make sense).  This term “dell-dot-F” then tells us the divergence of the F function.  If F were three directional flux content, this would be the amount of water leaving the point (x,y,z).  Notice that this is a scalar result!  So ∇F needs a scalar function F and gives a vector slope, while ∇·F uses a vector function F and gives a scalar result.


Icon for the Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License

Soil Hydrology and Biophysics by John Selker & Dani Or is licensed under a Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License, except where otherwise noted.

Share This Book