 
WIMOVAC Soil Processes
SoilPlant Water Relations
Models of water uptake by plant roots generally have one of two
purposes. Either they produce estimates of transpirational water loss for water budget
models or they provide estimates for predicting plant water stress {Campbell, 1991 #1958}.
Wimovac uses a simple computer algorithm to simulate multiple layer soil water flow, plant
water uptake and plantwater relations for a given soilplantatmosphere system. This is
based upon the analytical plant water uptake models derived by Campbell (1991) but
modified here to use wimovac’s own soil evaporation and canopy transpiration routines
in preference to the simpler formulations proposed by Campbell. Detailed reviews and
analyses of plant water uptake have been published {Molz, 1981 #1996; Hillel, 1980 #1995}
and form the background of this model module.
In the development of a theory to describe plant water uptake and loss,
it has been useful to employ electrical analogs of the system {Gardner, 1960 #1994; Cowan,
1965 #1992; Van den Honert, 1948 #1993}. Figure 33 shows the electrical analog used by
Cowan (1965) of a plant withdrawing water from a uniformly rooted soil layer (only
resistances in the liquid phase are shown).
Potentials and resistances in the
soilplant system. E is the transpiration rate, y_{L} is the leaf water potential,
y_{xL} is the xylem/leaf water potential_{ }, y_{xr} is the
xylem/root water potential, y_{r} root water potential, y_{s} is the soil
water potential, R_{L} is the leaf resistance, R_{x} is the xylem
resistance, R_{r} is the root resistance, R_{I} is the root interfacial
resistance and R_{s} is the soil resistance. Adapted from Campbell (1991)
Cowan (1965) assumes that transpiration rate (E) is regulated mainly by
vapour phase driving forces and resistances. Plant water potentials are therefore
primarily the result of an imposed rate of transpiration, from the leaf/canopy
transpiration module, and stomatal conductance rather than the converse situation in which
plant water potentials are the direct determinant of E. The major resistances in
Cowan’s model are outlined in Figure 33 and consist of resistances due to the soil (R_{s}),
soil interface (R_{I}), root endodermis (R_{r}), stem xylem (R_{x})
and leaf (R_{L}). For most species grown in well watered, medium textured soil,
the important resistances are due to the root endodermis (R_{r}) and the leaf (R_{L}).
In rough numbers a representative midday transpiration rate (E) of 2.3^{.}10^{4}
Kg m^{2} s^{1}, in a mature plant, produces leaf water potentials of
about 1500 kPa in exposed leaves and 1000 kPa in darkened leaves and from this
observation Cowan was able to parameterise Figure 33 such that when y_{s}=0,
y_{L}=1500 kPa^{ }and y_{xL}=1000
kPa. Cowan (1965) then assumed that the xylem and soil resistances were negligible and
could calculate R_{r} and R_{L} using Ohm’s law. Axial resistance in
the root, where important was assumed to be included in the R_{r} value. The root
resistance calculated in this manner is the total resistance of the root system. If the
length of root is known the resistance per unit length can be calculated and typical
values of around 2.5^{.}10^{10} m^{3} kg^{1} s^{1}
per unit length have observed {Bristow, 1984 #1997}.
Following the analysis of Cowan (1965) it is then possible to write a
differential equation for the uptake of water by a single root as:
where q_{w} is the flux of water (kg s^{1} ), k is the
soil hydraulic conductivity (kg s m^{3}), y is the
matrix potential (k Pa) and r is the radial distance from the root axis (m). The area for
water flow is A_{w}=2prl where l is the length of root.
Hydraulic conductivity can then be related to water potential using the term
where k_{s} is saturated conductivity (kg s m^{3}), y_{e} is the air entry potential (J kg^{1}) and n is
a constant that depends on soil texture and ranges typically from 2 to 3.5 {Campbell, 1974
#1998}. It is then possible to combine and integrate Equation 132 and Equation 133 from
the bulk soil at r_{s}, where water potential is y_{s},
to the root surface (r_{r}) where the potential is y_{r},
to give an expression (Equation 134) that may be related to measurable values.
where L is the rooting density and Dz is
rooting depth. Using this and the expression of Gardner (1960) {Gardner, 1960 #1994}for r_{s}
it is possible to calculate the soil resistance term as:

Equation 135 
Since B depends on root depth and root density, soil resistance in this
model is therefore a function of these variables. It will also vary with soil water
potential and to some degree with transpiration rate (E).
Few natural rooting systems would fit Cowan’s initial assumption
of a uniform rooting density and this model in itself would be of little use. However
Campbell (1991) has been able to extend this simple model to include variable rooting
density conditions in which the soil profile is assumed to be made up of a series of zones
or layers each of which has a characteristic root density, soil water potential and
associated resistances and it is this model that has been incorporated into wimovac.
In Campbell’s (1991) model water is assumed to move from the soil,
through the plant, to the evaporating surfaces in the substomatal cavities in response to
gradients in water potential. Only matric and gravitational potential are explicitly
accounted for in the derivation, and resistances to flow in the liquid phase are assumed
to be mainly in the soil, the root endodermis and the leaf. These resistances are assumed
constant and known during the period of a single simulation time step (15 minutes).
Interfacial and xylem resistance could have been included but have been assumed to be
negligible for the purpose of simplification in the model.
The water uptake from each soil zone in Campbell’s model is
described by the equations previously used by Cowan (1965) to describe a single layer
since the rooting density within a single zone can be assumed to be constant. The total
uptake of water from the soil is the sum of the uptake from each zone (Equation 136) where
y_{si} is the soil matric potential of zone i, g is the
gravitational constant (9.8 m s^{2}), z_{i} is the depth of layer i, y_{x} is the xylem water potential, R_{si} is the
soil resistance in zone i and R_{ri} is the root resistance calculated from
Equation 137 in which L_{i} is the rooting density of soil zone i. Equation 137
expresses the underlying assumption of the Campbell model that R_{ri} is directly
proportional to the total root resistance and inversely proportional to the fraction of
the root system that is in that layer. It is possible to then calculate the crown water
potential of the plant (y_{x}) using Equation 138 and
to use this to calculate the leaf water potential y_{L}
using Equation 139. y_{L }must then be solved
simultaneously with another equation that expresses the effect of y_{L}
on E through stomatal conductance (Equation 43, Equation 44) and this is achieved using an
iterative process.
Figure 34 shows typical diurnal output from the soilplant water
potential model for a an initially well watered clayloam soil during a 10 day dry down
period and the values generated for both leaf and soil layer water potentials are in
keeping with a number of field observations {Turner, 1983 #2000; Campbell, 1976 #1999}.
The diurnal patterns of leaf and soilwater potential variation are important because they
drive the stomatal closure model (Equation 43) which in turn strongly influences canopy
transpiration (E), assimilation (A) and energy budget. The short time step that the model
operates with (between 15 minutes and 2 hours) allows this model to express midday stress
on stomatal conductance with the consequent depression in midday transpiration and
assimilation typically observed in field trials {Campbell, 1976 #1999}.
Figure 34iii shows that the weighted mean soilwater potential for all soil layers is
strongly influenced by the densely rooted layers at the top of the soil profile. The fact
that mean water content (Figure 34iv) increases during the night means that the water
potential in these layers increases overnight as well. Campbell’s model does not
artificially restrict movement of water from the roots to the soil and assumes that that
resistance to water flow is the same for water flowing out of roots as it is for water
flowing in. This assumption is supported by work by Baker and van Bavel (1986){Baker, 1986
#2001}.
Soilplantwater relations. i)
Predicted decline in average leaf water potential during a 10 day dry down period
(British, midsummer) starting from a well watered condition. F_{canopy}=3. ii)
Increased stomatal closure resulting from lowered leaf water potential iii) Decline in
soil water potential of 5 uppermost soil layers (total soil depth represented by layers
1m). iv) Decline in volumetric water content of the 5 uppermost soil layers. All parameter
values as stated in Campbell (1976) unless otherwise stated in appendices I & II.
In addition to the relatively complex soilplant model of Campbell
(1976) wimovac also contains a much simpler model of volumetric soil water content
dynamics as proposed by Johnson (1993). {Johnson, 1993 #1919}
In his approach Johnson (1993) assumes a simple 1 layer soil system in
which there is no evaporation from the soil, runoff or drainage. The basis of the model is
to calculate the potential daily transpiration rate (the transpiration rate that occurs
when there is no resistance to water vapour transfer between the soil and the atmosphere)
using either the PenmanMonteith or the PriestlyTaylor formulations (Equation 146 and
Equation 147 respectively) and then to assume that the actual daily transpiration is
related to soil water content (q_{soil}).
For values of q_{soil}
greater than a particular critical value, q^{1}, (which
you can define), transpiration occurs at the potential rate, and for values of q_{soil} less than another
critical value, q^{2 }(the wilting point), there is no
transpiration. Intermediate values of q_{soil}
lead to transpiration rates linearly distributed between zero and the potential rate. Any
soil type can therefore be described in this model using three simple parameters, i). the
field capacity (q^{*}). ii). the critical value (q^{1}). iii). and the wilting point (q^{2}).
The default values for these parameters are included in wimovac in the form of a small
soil database that includes entries for clay, clayloam, sandy loam and user defined soil
types. The daily soil water balance is given by Equation 141 where the subscripts i and
i+1 refer to days i and i+1, r_{w} is the density of
water and d_{s} is the soil depth. Simplicity of parameterisation and design are
the main advantages to this approach but the direct, albeit empirical, manner in which the
rate of soil water loss is linked to soil water content provides a reasonable simulation
of the observed reduction in water loss rate associated with soil dry down evident in a
number of field observations (Campbell, 1976).
Equations from Campbell (1974) and Cowan (1965){Campbell, 1974 #1998; Cowan, 1965 #1992}.
Soil surface evaporation
The evaporation of water directly from the soil is often an important
component of effective water loss for the canopy, but is strongly dependent on soil
wetness and on plant cover {Jones, 1992 #1954}. Several reports indicate that, even with
wet soil, soil evaporation is only about 5% of the total when the leaf area index reaches
4 or more, for crops as diverse as wheat and Pinus radiata {Brun, 1972 #1894;
Denmead, 1969 #1893}{Denmead, 1969 #1893}{Brun, 1972 #1894}. When the leaf area index is 2
or less, however, wet soil evaporation can be as much as half the total. Ritchie
(1972){Ritchie, 1972 #1990} has suggested that after rain about 10mm of soil evaporation
can occur at a high rate (i.e. the soil surface resistance is close to zero). As the soil
dries, the surface resistance increases and soil evaporation decreases significantly. With
a dry soil surface, soil evaporation is relatively unimportant, and evaporation is solely
transpiration from plant leaves and this is approximately proportional to leaf area index
{Jones, 1992 #1954}.
Wimovac incorporates two alternative models of soil evaporation. i).
The first is commonly called the Priestly/Taylor approach and calculates a potential
evaporation rate by assuming that there is no resistance to water flux from the soil
surface. This is equivalent to the PenmanMonteith equation for potential transpiration if
the diffusion term is 26% of the radiation value. There appears to be little theoretical
basis for this simplification, but it is widely used by many modellers and is included
here for comparison. ii). The alternative soil evaporation model is based upon the
PenmanMonteith formulation which is similar in structure to that used for the calculation
of canopy transpiration described previously. This formulation is generally considered
more realistic and uses a soil boundary layer conductance term to give actual evaporation
rather than the potential evaporation produced by the Priestly/Taylor method.
Both of the soil evaporation models calculate the total incident solar
radiation on the soil surface and use this to calculate the net radiation balance (F _{N,soil}) of the soil by subtracting the longwave
radiation reemitted, from the soil, from the total incident radiation (Equation 145). For
the sunlit/shaded leaf canopy model, total incident solar radiation for the soil is
assumed to be equal to the leaf area weighted average light intensity incident on sunlit
and shaded leaves, and for the multiple layer canopy it is assumed to be equal to the leaf
area weighted average light intensity of sunlit and shaded leaves in the bottom canopy
layer. The radiation energy loss due to reemission from the soil is assumed to be given
by Stefans law reworked for an emitter and air at a similar temperature (Equation 143).
The soil boundary layer conductance g_{a,soil} is calculated using an empirical
expression (Equation 142) containing a term for the boundary layer thickness () where (u_{soil}) is the wind speed at the soil surface
and (S_{size}) is the average soil particle size and a diffusion component .
Both soil evaporation models show the expected reduction in soil
evaporation rate associated with larger canopy leaf areas but neither adequately takes
into account the increase in effective soil resistance that occurs during top layer soil
drying. Consequently both models may be prone to overestimate soil evaporation under dry
soil conditions. The PenmanMonteith boundary layer conductance (Equation 142) term could
be modified to include a soil resistance factor produced by the soil water potential
module but this is beyond the scope of current work and may represent an unnecessary
complication.
Equations from Johnson (1993).{Johnson, 1993 #1919}.
Heat flux & temperature
Models of soil heat flux and associated temperature estimates typically
make two major contributions to soil and crop models. Firstly they describes energy
partitioning at the soil surface which directly and indirectly affects plant growth and
development, and secondly they describe the soil temperature distribution which affects
many of the biological and physical processes associated with the soil.
Wimovac uses a one dimensional multiple layer model of soil temperature
and heat fluxes identical in structure to that used for the soil water potential
calculation module. There are four methods of generating soil temperature and heat flux
information included in wimovac and these may be selected by using the appropriate switch
in the soil microclimatetemperature section of the model parameter database. These
options reflect a series of trade off’s between ease of parameterisation and
mechanistic detail. i). The first option simply sets all soil layers to a given fixed
temperature which may be specified at the start of the simulation run in the model
parameter database. ii). The second option allows each soil layer to have a different
temperature fixed at the beginning of the simulation. iii). The third option sets all soil
layer temperatures equal to the ambient air temperature which changes dynamically with
time of day and day of year. There is little theoretical or experimental justification for
any of these options but they represent simple methods of specifying soil temperature and
are useful for debugging other complex model components which require an input of soil
properties. iv). The final option is by far the most mechanistically rich approach and
attempts to use the physically correct model proposed by Horton and Chung (1991).
Soil is represented in this model as a complex 3 phase system in which
there are solid particulates, a liquid phase water based flux and an associated vapour
flux arising from various physical and biological processes. Each type of soil represents
a potentially different balance of these components and a successful model of soil heat
flux and temperature must take these into account. The model proposed by Horton and Chung
(1991), and used in wimovac, assumes that {Horton, 1991 #1959} there are three mechanisms,
radiation, convection, and conduction responsible, simultaneously, for the transfer of
heat in soil. Radiative energy transfer includes incoming direct and diffuse shortwave
solar radiation, longwave sky radiation to the soil surface, and longwave radiation
reemitted outward from the soil surface. Convective energy transfer in a porous media
such as soil is associated with a net flux of fluids. A number of authors {Jackson, 1960
#2003; Wierenga, 1969 #2004} have studied the convective energy transfer in soils
associated with vapour and liquid fluxes in the soil and although convection may be a
major portion of heat flow during periods of large moisture flux, such as during rainfall
or irrigation, the mechanism most responsible for subsurface heat transfer is conduction
{Horton, 1991 #1959}.
Philip and De Vries (1957) {Philip, 1957 #2008}presented a theory to
describe coupled water and heat flow in soil. Van Bavel and Hillel (1975, 1976){Van Bavel,
1976 #2010}{Van Bavel, 1975 #2009}, Milly (1982) {Milly, 1982 #2011} and Bristow et al.,
(1986) {Bristow, 1986 #2012}expanded and/or used this theory to calculate heat and
water flows in soil and were able to derive a description (Equation 148) for conductive
heat transfer for a vertical one dimensional soil system.
where C_{soil} is the volumetric heat capacity (J ° C^{1} m^{3}), T_{a} is the
temperature (° C), t is the time(s), z is the
depth (positive downward, m), and l_{soil} is the
thermal conductivity (W m^{1} ° C^{1}). In
the presence of a temperature gradient in a moist soil, heat transfer takes place by
convection, in addition to conductance, but by making l_{soil}
the apparent rather than the real thermal conductivity Horton and Chung (1991) avoid
having to include this term explicitly. In a moist soil l_{soil
}depends on the soil water content which can vary with both depth and time.
The movement of water in the one dimensional soil profile is described
by Equation 149.
where F is the specific water capacity, ,
and where q_{soil} is the volumetric soil water
content, h_{soil} is the soil water pressure head, K_{soil} is the
hydraulic conductivity, S is the source or sink term and z_{soil} is the vertical
distance downward from the soil surface.
Horton and Chung (1991) showed that Equation 148 and Equation 149 can
then be solved simultaneously provided that the initial and boundary conditions are known.
The initial conditions in the wimovac version of this model are
specified in the model parameter database and consist of initial soil layer temperature
information and various soil physical properties (see appendix II). Boundary conditions
are calculated on the basis of an energy balance method at the soil surface which is used
to determine thermal and hydraulic upper boundary properties. Soil evaporation is
calculated in the soil surface evaporation module (Equation 142 to Equation 147) and is
used in conjunction with incident and reemitted solar radiation in order to calculate the
latent and sensible heat fluxes at the surface. Soil surface evaporative flux and soil
surface temperature are not independent quantities and so are estimated in the model in an
iterative fashion. The humidity at the soil surface is calculated from Equation 150, where
HO_{soil} is the saturated humidity at the soil surface (kg m^{3}) and h_{soil}
is the soil water pressure head.
The saturated humidity, HO_{soil}, at the soil surface is
calculated using Equation 151 derived by Murray (1967) {Murray, 1967 #1879} and the soil
heat flux density, G_{soil}, at the soil surface may then be described by Equation
152 where l_{soil} is the thermal conductivity of the
soil surface (W m^{1} °C^{1}). Net radiation, soil evaporation, soil
surface humidity and G_{soil} are all functions of the unknown soil surface
temperature and before the surface temperature can be calculated an approximation for G_{soil}
must be used (Equation 153) , this is provided by Chung and Horton (1987){Chung, 1987
#2013} in which T_{s} is the surface temperature for the present time step, T_{1}
is the surface temperature from the previous time step, T_{2} is the soil
temperature from the previous time step for the node at vertical position 2, Dz is the vertical spatial increment, C is the volumetric heat
capacity for the soil surface layer and Dt is the time step
increment.
Equation 153 approximates the soil heat flux density by adding a term
that estimates soil heat flux at a depth of and a term that estimates the change in heat stored in the soil above . The value of
T_{s} is solved numerically for each time step using a bisector rootfinding
algorithm {James, 1977 #2014}. The predicted soil surface temperature, T_{s} from
the energy balance partitioning is used as the upperboundary condition. The model
algorithm then uses an implicit, finite difference method (Equation 154) to solve the soil
heat and water flow equations.
Soil temperature profiles. i).
Predicted soil temperature profile for J_{day}=190 at 0400hr and 1600hr indicating
the change in temperature in the upper layers. Note that the upper levels are more
directly connected to air temperature compared to lower layers which are relatively more
stable. ii).Predicted soil temperature surface for J_{day}=100 to J_{day}=101.
Total depth of soil profile=0.5m, F_{canopy}=3. Parameters values as stated in
Horton and Chung (1991) unless otherwise specified in appendices I & II.
Soil parameters required by the model as inputs include soil surface
emissivity, initial soil surface albedo, soil thermal conductivity, soil volumetric heat
capacity, soil water capacity, and hydraulic conductivity, all as functions of water
content, roughness length, and the soilwater characteristic curve. A few general inputs
are also required to run the module and these are obtained from the appropriate routines
in the macroclimate module and include the time of solar noon, day length, length of
simulation, D z, and D l and various
weather inputs such as rainfall. Soil water characteristics, hydraulic conductivity
and specific water capacity are described by empirical equations presented by Van
Genuchten (1980) as given by Equation 155 in which K_{s} is the saturated
hydraulic conductivity, h is the absolute value of the pressure head and a _{1} and n are non linear regression parameters describing
the shape of the soil water characteristic curve. Soil surface emissivity, e , (Equation 157) follows that used by Van Bavel and Hillel
(1976).{Van Bavel, 1975 #2009} The soil volumetric heat capacity (J °C^{1} m^{3})
is determined following De Vries (1963) {De Vries, 1963 #2015}and is given by Equation
158, where q _{s} is the input saturated volumetric
water content of a specified soil. Soil thermal conductivity is assumed to be related to
water content by the empirical expression (Equation 159). Soil surface albedo (a_{soil})
may be obtained as a fixed value from the model parameter database or determined according
to Van Bavel and Hillel (1976) as indicated in Equation 160.
Equations from Horton and Chung (1991), {Horton, 1991 #1959}
